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Chapter  1 
INTRODUCTION 

Increased  competition,  coupled  with  the  recessive  climate  of 
today's  marketplace,  require  effective  control  and  maintenance  of 
inventory  systems.   The  proper  balance  between  inventory  investment  and 
profits  can  make  the  difference  between  the  success  or  failure  of  a 
business  (1).   The  fundamental  reason  for  maintaining  items  in  inventory 
is  that  it  is  physically  impossible  or  economically  unsound  to  have  an 
item  arrive  at  the  precise  instant  a  demand  for  that  item  occurs. 
Businesses  thus  typically  carry  in  inventory  raw  materials,  finished 
goods,  and  repair  or  replacement  components. 

The  manufacturing  industries  realized  the  need  for  inventories  in 
the  late  1800s  because  finished  products  were  being  produced  in  lots, 
production  set-up  costs  were  high,  and  storage  of  the  finished  goods  in 
factory  warehouses  was  necessary.   Prior  to  1940  the  inventory  systems 
were  treated  as  deterministic  systems  since  component  needs  were  assumed 
to  be  known  quantities  (2) .   The  simple  lot-size  formula  was  derived 
first  by  Harris  and  later  by  Wilson,  who  used  it  to  analyze  inventory 
systems  (3).   During  the  early  1940s,  the  "Christmas  tree  model"  or  the 
"newsboy  problem"  was  the  first  attempt  to  model  stochastically  a  demand 
for  a  product  (4).   After  1945,  the  newly  developed  fields  of  management 
science  and  operations  research  also  investigated  the  stochastic  nature 
of  inventory  problems. 

Once  of  the  early  analytical  models  developed  to  describe  inventory 
systems  was  queueing  theory.   A  queue  or  waiting  line  is  formed  by 
randomly  generated  demands  at  a  service  facility.   Queueing  theory  can 
be  used  to  predict  the  steady-state  probability  of  a  system  being  in  a 


particular  state  at  a  specified  time.   The  pioneering  work  using 
queueing  theory  was  performed  by  Erlang  during  1901  to  1920  (5),  while 
the  stochastic  treatment  of  queueing  theory  was  studied  by  Kendall 
during  1951  to  1953  (6,7).   Since  that  time  numerous  studies  have  been 
conducted  using  queueing  theory  to  model  inventory  systems.   Implicit 
assumptions  of  queueing  theory  are  an  infinite  calling  source  and  the 
independence  of  different  states  of  the  system  (8,9). 

There  is  considerable  incentive  to  evaluate  the  economic  feasibility 
of  spare  pool  inventory  systems  for  high  cost  components.   A  stochastic 
process  is  used  to  determine  the  demand  on  the  pool  by  combining  two 
frequently  researched  areas  of  inventory  systems;  namely,  the 
stockpiling  of  low  cost  components  and  the  repair  or  replacement  of  high 
cost  components.   Low  cost  components  are  usually  stockpiled  if  their 
failure  could  disrupt  the  system  operation.   Upon  failure  of  a  low  cost 
component,  the  failed  component  is  replaced  with  a  new  component  and 
then  discarded.   The  loss  of  production  is  insignificant  if  the 
installation  time  is  short.   By  contrast,  high  cost  components  are  not 
usually  held  in  inventory  and  must  be  repaired  or  replaced  upon  failure. 
The  loss  of  production  can  be  minimized  by  shifting  manpower  and  other 
resources  to  another  production  line  during  the  repair/replacement 
and/or  installation  time. 

A  different  situation  arises  for  certain  component  failures  in  an 
electrical  power  plant  where  the  only  product  is  electricity.   During 
component  failures  that  result  in  shutdown  or  derating  of  the  plant, 
manpower  and  other  resources  cannot  be  shifted  to  an  alternate 
production  line.   Of  particular  interest  are  components  with  low 
probabilities  of  failure  and  large  capital  costs  (e.g.,  typically 


hundreds  of  dollars).   Main  power  transformers,  coolant  water  pumps, 
turbine  rotors,  and  nuclear  safety  equipment  failures  are  examples  of 
components  that  can  cause  a  significant  loss  of  output  capacity. 
Compounding  the  problem  is  the  extended  length  of  time  a  plant  could  be 
shut  down  or  derated  because  of  a  long  repair/replacement  time  for  these 
components . 

A  management  scheme,  known  as  Pooled  Inventory  Management,  was 
developed  by  the  General  Electric  Company  (GE)  to  help  prevent  extended 
outages  at  nuclear  power  plants  (10).   In  this  program,  components  with 
small  failure  probabilities  and  long  repair/replacement  times,  and  whose 
failure  would  result  in  a  significant  reduction  in  plant  capacity  are 
placed  in  a  spare  pool.   The  costs  of  the  spare  pool  operation  are 
shared  among  the  members  of  the  pool  thereby  reducing  the  added  cost  of 
the  spare  purchase  while  providing  the  advantages  of  a  spare  component 
to  each  pool  member. 

The  economic  evaluations  of  such  a  spare-component  pool  is 
generally  quite  complex  as  a  result  of  differing  economic  procedures 
used  by  pool  members,  differences  in  the  capacity  of  members'  plant 
components,  and  the  stochastic  nature  of  the  failure  process.   A  1974 
study  conducted  by  GE  proposed  a  method  to  evaluate  a  spare-component 
pool  for  generator  step-up  transformers  at  electrical  power  stations 
(11).   While  this  study  was  somewhat  simplistic  (e.g.,  infinite  pool 
lifetime,  no  escalation  of  costs  and  all  costs  are  capitalized),  it  did 
demonstrate  the  economic  attractiveness  of  sharing  a  spare  step-up 
transformer  among  several  utilities.   In  1985,  Shultis  (12)  performed  a 
more  detailed  analysis  of  the  spare-pool  problem,  in  which  two  component 
management  plans,  designated  by  Plan  A  and  Plan  B,  were  evaluated. 


Under  Plan  A,  each  plant  repairs  or  replaces  its  failed  component  when  a 
failure  occurs.   Under  Plan  B,  a  spare  component  is  purchased  and  placed 
in  storage  prior  to  failure.   If  at  the  time  of  a  component  failure  at  a 
plant  the  spare  is  available,  the  failed  component  is  removed  and 
replaced  by  the  spare.   The  spare  is  used  as  a  temporary  substitute 
until  the  failed  component  is  repaired  or  replaced. 

Shultis  estimated  the  probability  of  spare  availability  using 
binomial  theory  to  calculate  the  fraction  of  each  year  the  spare  is 
available.   As  a  consequence,  the  spare  is  used  only  if  it  is 
immediately  available  at  the  time  of  a  plant  failure  and  the 
repair/replacement  time  is  less  than  one  year.   Because  the  components 
considered  in  this  study  have  very  small  failure  rates,  the  probability 
of  multiple  failures  of  a  single  component  was  considered  negligible. 
Similarly,  the  probability  of  the  spare  failing  during  the  temporary 
installation  time  was  neglected.   In  addition,  the  economic  costs  are 
estimated  on  a  yearly  discounted  basis  with  the  assumption  that  all 
costs  are  paid  at  the  end  of  each  yearly  interval. 

The  purpose  of  the  present  study  was  to  refine  some  of  the 
probabilistic  and  economic  models  used  by  Shultis  and  to  assess  the 
importance  of  such  refinements.   Specific  objectives  were  to: 

(a)  develop  a  general  Monte  Carlo  method  to  describe  the  spare 
component  availability  and  fraction  downtime  (i.e.,  fraction 
of  the  repair /replacement  time  the  plant  is  shutdown  or 
derated) , 

(b)  develop  a  model  for  continuous-time  cost  analysis  in  which 
costs  are  evaluated  at  the  time  they  are  incurred  and  not  at 
the  end  of  the  year  as  was  done  in  the  earlier  study, 

(c)  compare  how  different  failure-model  assumptions  affect  the 
subsequent  economic  analysis  of  the  spare-component  pool,  and 

(d)  compare  results  obtained  via  the  continuous-time  cost  model 
with  the  year-end  cost  model. 


As  in  the  previous  study,  total  costs  are  estimated  for  two 
component  plans,  designated  by  Plan  A  and  Plan  B.   Under  Plan  A,  each 
plant  repairs  or  replaces  the  failed  component  when  a  failure  occurs. 
The  plants  are  assumed  to  be  operating  independently  in  the  sense  that  a 
failure  in  one  plant  has  no  effect  on  another  plant.   In  the  case  of  a 
component  failure,  the  plant  will  remain  shutdown  or  derated  for  a  fixed 
period  of  time.  There  is  no  provision  for  an  operating  plant 
with  a  failed  component  to  borrow  a  substitute  component  from  another 
plant. 

Under  Plan  B,  a  single  spare  component  is  purchased  and  placed  in 
storage  prior  to  component  failure.   If  at  the  time  of  a  component 
failure  at  a  plant  the  spare  is  available,  the  failed  component  is 
removed  and  replaced  by  the  spare  and  used  as  a  temporary  substitute 
until  the  failed  component  is  repaired  or  replaced.   The  plant  with  the 
failed  component  has  the  option  to  use  the  spare  should  it  become 
available  prior  to  the  arrival  of  the  repaired/replacement  component. 
If  the  spare  is  installed  during  the  repair/replacement  time  the  plant 
avoids  a  portion  of  the  costs  that  would  have  to  be  paid  under  Plan  A. 
However,  the  plant  must  also  pay  for  a  portion  of  the  spare-pool  costs 
in  addition  to  the  expenses  incurred  during  removal  and  installation  of 
the  failed  component  and  spare.  During  the  time  the  spare  is  not 
installed,  the  failed  plant  incurs  costs  as  if  operated  under  Plan  A. 

This  study  considers  only  the  case  when  the  repair/replacement  time 
is  a  fixed  period  of  time.   In  addition,  the  repair/replacement  time  is 
assumed  to  be  the  same  for  each  plant  with  a  similar  component 
regardless  of  the  component's  capacity  or  operating  history. 


For  the  Plan  B  analysis,  the  fraction  of  the  repair/replacement 
time  the  plant  is  shut  down  or  derated  due  to  the  failure  of  the 
component  under  consideration  was  designated  as  the  "fraction  downtime". 
Because  of  the  addition  of  a  spare  to  the  system  and  a  finite  number  of 
plants  in  the  pool,  standard  queueing  theory  models  cannot  be  used  to 
estimate  the  fraction  downtime.   An  exact  analytical  solution  for  the 
fraction  downtime  that  takes  into  account  all  possible  combinations  of 
pool  lifetimes,  removal  and  installation  times,  component 
repair/replacement  times,  etc.  appears  to  be,  if  even  possible,  very 
difficult.   However,  because  of  the  uncertainties  in  the  component 
failure  rates  and  the  simplified  engineering  economic  analysis  used,  the 
complexity  involved  in  a  detailed  theoretical  approach  may  not  be 
warranted . 

An  approximate  solution  for  the  fraction  downtime  can  be  obtained 
using  renewal  theory  techniques  (13).   In  a  system  described  by  standard 
renewal  theory,  a  plant  can  be  in  only  one  of  two  states;  that  is,  a 
plant  operating  with  an  original  or  repaired/replacement  component  is 
said  to  be  in  the  "on"  state.   A  plant  that  is  down  or  operating  with 
the  spare  installed  is  in  the  "off"  state.   The  spare  is  available  at 
the  time  of  a  plant  failure  only  if  all  plants  in  the  system  are  in  the 
"on"  state  excluding  the  "just  failed"  plant.  Standard  renewal  theory 
is,  however,  limited  to  estimating  the  fraction  downtime  for  the  case 
when  the  spare  is  used  only  by  the  "just  failed"  plant  if  immediately 
available  at  the  time  of  failure  and  the  spare  failure  rate  is  zero. 

Another  approximate  method  to  estimate  the  expected  fraction 
downtime  can  be  obtained  by  simulating  the  component-failure  problem 
using  a  Monte  Carlo  analysis  (14,15).   An  approach  similar  to  that  used 


by  Widawsky  (14)  and  Emshoff  and  Sisson  (15)  was  used  where  the  variable 
time-increment  method  was  used.  The  random  failure  times  were  simulated 
by  generating  a  sequence  of  random  numbers  from  the  appropriate 
probability  distributions  characteristic  of  the  desired  failure  model. 
The  simulated  failure  times  for  each  pool  history  were  then  analyzed  to 
determine  when  the  spare  would  be  available  for  use  as  well  as  the 
fraction  of  the  repair/replacement  time  for  which  the  spare  could  be 
installed  at  a  failed  plant.   A  large  number  of  pool  histories  were 
simulated  and  the  results  of  each  history  were  then  averaged  to  estimate 
the  probability  of  spare  availability  and  the  fraction  downtime  given  an 
initial  set  of  conditions. 

A  computer  code,  SIMULATION,  was  written  based  on  the  Monte  Carlo 
simulation  technique  for  three  different  failure  models.   Model  1  had 
the  following  properties:   (i)  the  spare  has  zero  failure  probability, 
(ii)  the  components  are  installed  and  removed  instantaneously,  and  (iii) 
the  spare  is  used  when  it  becomes  available.   Model  2  was  identical  to 
Model  1  with  the  exception  that  the  spare  is  used  only  if  it  is 
available  at  the  time  of  a  component  failure.   Model  3  was  also  similar 
to  Model  1  with  the  exception  that  the  spare  has  a  failure  rate  equal  to 
the  failure  rate  of  the  components  installed  in  the  plants. 

To  verify  the  computer  simulation  the  results  of  Model  1  were 
compared  to  the  results  predicted  by  renewal  theory.   Renewal  theory  can 
determine  the  probability  that  a  system  is  in  either  an  "on"  or  "off" 
condition.   For  the  case  of  a  spare  with  zero  failure  probability  and  a 
fixed  repair/replacement  time  (Model  1)  the  probability  of  spare 
availability  at  a  given  time  to  that  estimated  by  the  computer 
simulation  should  agree  with  the  probability  that  all  plants  are 
operating  at  time  t  as  estimated  by  renewal  theory. 
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As  a  second  phase  of  evaluating  the  effectiveness  of  a 
spare-component  pool,  the  failure  models  were  then  combined  with 
economic  models  to  determine  the  present  worth  costs  of  the  two 
component-management  plans.  A  simple  economic  analysis  approach  was 
used  that  allowed  for  yearly  varying  escalation  and  discount  rates  (16). 
Unlike  similar  analyses,  the  economic  model  formulated  in  this  research 
project  estimated  present  worth  costs  based  on  when  the  costs  were 
incurred  instead  of  the  typical  discrete  year-end  method. 

The  economic  model  used  in  this  study  divided  the  costs  into  four 
major  groups.   The  failure-dependent  variable  costs  included  the  cost  of 
purchasing  replacement  energy  during  a  failure.   The  variable  costs 
during  plant  operation  included  the  operation  and  maintenance  costs 
associated  with  an  operating  component.   The  failure-dependent  fixed 
costs  included  the  repair/replacement  component  costs,  the  failed 
component  salvage  value,  and  the  cost  of  temporarily  installing  the 
spare.   The  annual  fixed  costs  included  reserve  capacity  costs, 
spare-pool  maintenance  costs,  spare  salvage  value,  and  the  used 
component  salvage  value. 

A  computer  code,  KSUSPARE,  was  written  to  calculate  the  estimated 
Plan  A  and  Plan  B  present  worth  costs.   The  analytical  methods  used  in 
KSUSPARE  were  then  compared  to  the  analytical  methods  used  in  PC-SPARE. 
The  fraction  downtime  estimated  by  SIMULATION  for  all  three  simulation 
models  and  renewal  theory  were  used  in  KSUSPARE  to  estimate  Plan  A  and 
Plan  B  present  worth  costs.   Finally,  an  example  problem  compared  the 
numerical  differences  between  PC-SPARE  (12)  and  KSUSPARE  for  three 
variations  of  initial  conditions. 


The  present  study  is  divided  into  three  major  sections.  Chapters 
II  and  III  present  the  development  of  the  failure  and  economic  models, 
respectively,  while  in  Chapter  IV  the  continuous-time  and  year-end  cost 
models  are  compared.   Chapter  V  summarizes  the  study  and  suggests  areas 
for  further  study. 
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Chapter  2 
FAILURE  MODELS 

To  assess  the  economic  feasibility  of  creating  a  spare-component 
pool,  it  is  first  necessary  to  estimate  the  expected  or  average  downtime 
for  a  given  plant  over  the  spare-pool  lifetime.   The  downtime  estimate 
includes  the  interdependence  of  the  possible  failure  times  for  all 
plants  in  the  spare  pool  and  availability  of  the  spare  in  the 
repair/replacement  time  interval.   In  this  study  the  repair/replacement 
time  for  the  failed  component  is  a  fixed  time.   Simplifying  assumptions 
are  that  similar  plant  components  have  equal  failure  rates  and  equal 
repair/replacement  times. 

An  exact  analytical  solution  for  the  average  downtime  that  accounts 
for  all  possible  combinations  of  pool  lifetimes,  component  removal  and 
installation  times,  component  repair/replacement  times,  etc.  appears  to 
be  very  difficult.  For  a  negligible  spare  component  failure  rate  and 
spare  usage  assumptions,  the  average  downtime  can  be  estimated  by 
renewal  theory.   The  average  downtime  for  the  more  general  case 
conditions  must  be  estimated  using  Monte  Carlo  techniques. 

2.1  Plant  Operation  as  a  Renewal  Process 

Consider  a  single  plant  with  a  component  failure  rate  A  when  the 
original  component  operates.   The  plant  is  either  "on"  (i.e.,  operating) 
or  "off"  (i.e.,  shutdown).   A  system  that  can  be  described  as  being  in 
one  of  two  states,  "on"  or  "off",  is  a  renewal  process. 

A  renewal  process  is  a  series  of  points  on  the  interval  (0, °$   where 

the  times  between  failures  are  independent  and  identically  distributed 

random  variables  (17).   Let  T   be  a  random  variable  defined  as  the  time 

on 
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before  failure  of  the  plant  component.   The  cumulative  distribution 

function  (CDF)  for  the  random  variable  T   is  given  by 

on    ° 


F   (t)  =    f   (x)dx,  (2.1) 

on      I    on 


where  f   (x)  is  the  associated  probability  density  function.   Let  T  ,_ 

be  the  random  variable  defined  as  the  repair/replacement  time.  The 

plant,  therefore,  operates  for  a  random  time  T   and  is  then  down  for  a 

random  time  T   .  Let  T  be  the  time  the  plant  begins  operating  after  a 

repair  (i.e.,  where  T  =  T   +  T  „,,  and  let  f(x)  be  its  density 
on    off 

function.   The  function  g(t)  is  defined  as  the  probability  the  system  is 
operating  at  time  t  given  by 

CO 

g(t)  =    Prob{plant  operating  at  t  |  T  =  x}  f(x)  dx.     (2.2) 
0 

For  the  case  x  <  t  the  probability  the  plant  is  operating  at  time  t 

is  equal  to  the  probability  the  plant  was  operating  at  time  t-x  (see 

Fig.  2.1).   For  the  case  x  >  t  the  probability  that  the  plant  is 

operating  at  time  t  is  equal  to  the  probability  that  t  is  in  the 

interval  [0,T   ].   Thus, 
on 

8 (t-x)  ,   x<t 

Prob{plant  operating  at  t  |  T=x)  =  I  (2.3) 

Prob{t  <  T    |  T=x}   ,   x>t. 


Substitution  of  this  result  into  Eq.  (2.2)  yields 
t  » 

g(t)  =    g(t-x)  f(x)  dx  +    Probk  <  T__  |        :(x)  dx.   (2.4) 


on 

0  t 


The  second  integral  in  the  above  equation  is  the  mathematical  expression 

for  Prob  {t  <  T   and  t  <  T}.   Because  T   <  T, 
on  on  — 
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Prob{t<T   |  T=x)  f(x)dx  -  Prob{t  <  T   and  t  <  T}  =  Prob{t  <  T  }. 
I  on  on  on 

t 

Thus,  the  results  from  the  above  expression  can  be  combined  with 

Eq.  (2.1)  to  give, 

t 

J   ProbCt  <  TQn  |  T-x}  f(x)dx  =  1  -  |  fQn(x)dx.        (2.5) 
t 

The  above  expression  can  now  be  combined  with  Eq.  (2.4)  to  give 
t  t 

g(t)  -  1  -  I  fon(x)dx  +  J  g(t-x)  f(x)  dx  .  (2.6) 

I   -co  J  0 

Equation  (2.6)  is  a  renewal  equation  for  which  expressions  for  f   (t) 

on 

and  g(t-x)  must  be  found. 


-T 
off 


r 


x=0  T        x=T      t 

on 

FIG.  2.1.   Graphical  representation  of  single  plant  operati 


on. 


Because  the  component  failure  rate  is  constant,  the  time  intervals 
between  each  component  failure  are  independently  and  identically 
distributed  according  to  the  negative  exponential  distribution  (17). 
The  value  of  x  is  defined  only  for  positive  values  of  time,  thus 

1  -  J   fon(x)dx  -  1  -  J   fQn(x)dx  =  1  -  J  X   e"Xx  dx  =  e"Xt.  (2.7) 
-"0  0 

The  Prob(T  <  t}  =  0  for  t  <  T    because  the  plant  must  have 

completed  a  full  operating  cycle  (i.e.,  T  +T  ,,)  before  operating 

on  off         r  ° 

again.   The  random  variable  TQn  can  equal  zero  but  the  minimum  value  of 
Toff  is  the  fixed  repair/replacement  time.  Thus,  the  CDF  for  T  is 
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(2.8) 


F(T)    =  Prob{T   <  t}  =  Prob{T        <   t-T  „} 
-  on  -         off 

J   1   -exp[-X(t-Toff)]      ,      t    >Toff    . 

The  density  function  f(T)  found  by  differentiation  of  F(T)  is 

0  .   t<   T 

f(T)  =  \  °"  (2.9) 

X  exp[-X(t-T   )]    ,    t  >   T    . 
ott  orr 


The  above  result  for  the  density  function  can  now  be  substituted  into 
Eq.  (2.7)  to  give 

.    t    <  T 

'        -     off 

-xt  f*         -A*  (2-10) 

+  X  exp(X  T        )  e   AX  g(t-x)dx    ,    T   ,,   <   t   <  ». 


g(t)    = 


off    L  *  5V"   *'—    '      off 

off 


2.2  Relationship  of  Renewal  Process  to  Spare-Component  Problem 

Consider  the  case  of  one  plant  in  the  spare  pool.   The  failure  rate 
of  a  one  plant  system  is  X^  when  the  original  component  operates.   If 
the  spare  component  is  perfect,  i.e.,  has  zero  failure  rate,  the  system 
failure  rate  is  zero  when  the  plant  is  shutdown  or  operating  with  the 
spare  installed.   For  these  specific  conditions  the  system  is  either 
"on"  (i.e.,  operating  without  the  spare  component)  or  "off"  (i.e.,  plant 
is  shutdown  or  operating  with  its  original  component) .   A  system  that 
can  be  described  as  either  "on"  or  "off"  is  a  renewal  process.   Thus, 
T  fc  equals  the  plant  component  repair/  replacement  time  during  which 
the  plant  is  shutdown  (or  operates  with  the  perfect  spare  as  a  temporary 
substitute).   Hence,  the  probability  that  the  plant  is  operating  in  its 


original  state  is  given  by  Eq.  (2.10)  with  T    =  T    . 
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A  substitution 


of  variables  can  be  made  where  x  =  t  -  u  and  dx  =  -  du.   Equation  (2.10) 
then  becomes 


g(t) 


-  * 


,  0  <t  <T 


repl 


t-T 


(2.11) 


♦    Xexp(ATrepl) 


repl 


exp[-x(t-u)]g(u)du,T     <t  < 


0 


2.3  Solution  Techniques  for  the  Renewal  Equation 

The  evaluation  of  Eq.  (2.10)  can  be  accomplished  by  using  three 
methods:   analytical,  numerical,  and  Laplace  transforms. 

2.3.1  Analytical  Solution  to  the  Renewal  Equation  for  the  First  Three 
Intervals 

The  renewal  equation  of  Eq.  (2.11)  can  be  cast  into  a  somewhat 

simpler  form  by  introducing  the  variables  g=  \T    ,  T=  Xt,  and 

t'  =  A  u.   Equation  (2.11)  then  becomes 


g(T)  = 


e"T+ee  j 


,    o  <  T  <  e 


e"(T"T,)  g(i')  dT'   ,   B  <  t  <  - 


(2.12) 


Let  g  (t)  be  the  solution  to  Eq.  (2.12),  where  g.(x)  =  gd)  for  (i-1) 
S  <  T  <  i  8.   The  solution  to  Eq.  (2.12)  for  the  first  three  T  intervals 
is  listed  below. 

First  interval  (0  <  t  <  1) 


g^t)  =  e 


(2.13) 


Second  interval  (1  <  X  <  2) 
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g2(T)  =  8i(t)  +  eB  I   e"(T"T,)  gjd'tdt', 


which,  upon  evalution  gives 

g2(t)  =  e"T  [l+e6(t-6)] 


Third  interval  (2  <  t  <  3) 


(2.14) 


T-B 
g3(T)  =  gl(T)  +  eP  j  I  e"(T_T')g1(T')dT'  +  f   e"(T"T,)g2(T')d(T') 


which  evaluates  to 


g3(t)  =  e   <1  +  e"  |t-S  ♦  eM  [(t/2-6)  (t-2)  ]  |  \  .  (2.15) 


The  above  integral  equation  demonstrates  the  dependence  of  g(i)  upon 
g,_.(t).   At  each  evaluation  time  T,  the  value  for  g(i)  depends  upon  the 
previous  values  of  g(T)  over  the  interval  (0,(i-l)t  ).   Thus,  the 
expressions  for  g.(x)  get  more  complex  as  i  increases.   An  alternative 
for  evaluating  g  (t)  analytically,  especially  for  large  t,  is  to  solve 
Eq.  (2. 12)  numerically. 

2.3.2  Numerical  Solution  to  the  Renewal  Equation 

The  method  of  numerical  integration  chosen  to  evalute  the  renewal 
equation  was  the  trapezoid  rule.   Consider  an  integrable  function  f(x) 
that  is  to  be  evaluated  over  the  interval  a  <^  x  <  b.   The  integral  I  is 
defined  as 


I  =    f(x)  dx. 


(2.16) 
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The  interval  [a,b]  is  then  divided  into  N  equal  subintervals  of  width 
&,  where 

*-TT-  (2-17) 


The  area  under  the  curve  in  each  subinterval  bounded  by  x    and  x  can 

be  approximated  by 

x.         „      . 

f J         f1-l  +  f1 

f(x)dx  s-J-t- J-  (Ac)  (2.18) 

where  f    and  f  denote  the  function  f(x)  evaluated  at  x.   and  x 
respectively.  Thus,  Eq.  (2.16)  can  be  approximated  by  (18) 

j  f(x)dx  *f-   (fj  ♦  fN+1  ♦  2  I  fj)  ,  (2.19) 


J-2 


where  fj  =  f(a)  and  f    -  f(b). 


The  trapezoid  rule  can  be  applied  to  the  renewal  equation. 
Starting  at  time  zero  each  repair/replacement  interval  8  is  split  into  N 
equal  subintervals  (see  Fig.  2.2). 


t-0  t=8  T»2B 

H+++H 1 K—time 

i-1    2   3  4  5    ...  N+l  2N+1 

FIG.  2.2.  The  time  line  broken  into  interval  lengths 
of  8  with  N  subintervals. 


From  Figure  2.2  it  can  be  seen  T  =  i  At  for  i  =  0,  1  ...  °°.   For  the 
case  of  i  <  N 
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-T.    -IAt, 


(2.20) 


while  for  i  >  N 


8(T1)  =  8i  =  e    i1 


eT  g(i')dt' 


(2.21) 


Applying  the  trapezoid  approximation  to  Eq.  (2.21)  gives 


8<  " 


1  ♦  e  y.  8l  ♦  e 


V-6 


1-8   Xt, 


H-N+l 


+  2 


j"2 


,  i  >  N.   (2.22) 


Because  the  analytical  solution  to  the  renewal  equation  is  known  for  X  < 
6,  the  first  intervals  of  i  up  to  N  are  calculated  using  Eq.  (2.20). 
The  remaining  intervals  are  calculated  using  the  results  from  i  up  to  N 
and  Eq.  (2.22). 

A  computer  code  was  written  to  evaluate  Eq.  (2.22).   The  numerical 
integration  was  performed  for  a  maximum  time  such  that  the  result  for 
g(t  >  4)  was  within  0.01%  of  the  asymptotic  value  (see  Section  2.3.3). 
The  result  from  the  numerical  integration  was  then  compared  to  the 
analytical  results  for  g(r)  derived  in  Eq.  (2.13),  (2.14),  and  (2.15) 
for  the  first  three  intervals.   The  results  indicate  excellent  agreement 
between  the  analytical  and  the  numerical  approximation  for  a  very  small 
number  of  subintervals.   It  is  this  numerical  procedure  that  was  used  in 
the  subsequent  simulation  of  plant  failures  outlined  in  Section  2.5. 

The  trapezoid  solution  technique  is  exact  when  compared  to  the 
analytical  solution  over  the  first  three  intervals.   The  g  (f)  functions 
being  integrated  over  these  intervals  are  equal  to  zero,  a  constant,  and 
a  linear  term,  respectively,  for  which  a  straight  line  representation  is 
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exact.  Thus,  a  fine  quadrature  mesh  is  necessary  only  for  intervals  at 
time  t  >  3. 

The  results  of  an  example  problem  for  a  single  plant  are 
represented  graphically  in  Figure  2.3.  A  study  of  this  figure  indicates 
a  negative  exponential  behavior  in  the  first  repair/replacement  time 
interval,  followed  by  oscillations  in  the  g(t)  value.   After  several 
repair/replacement  time  intervals  the  value  for  g(t)  stabilizes  to  an 
asymptotic  value.   Thus,  for  long  lifetimes,  the  approximation  can  be 
made  that  the  probability  a  plant  is  operating  at  time  t  equals  the 
asymptotic  solution  to  the  renewal  equation. 

2.3.3  Solving  the  Renewal  Equation  Using  Laplace  Transforms 

The  renewal  equation  for  our  one  component  system  can  be  rewritten 
as 

'  -Xt 
e 

g(t)    ■   '  t 

e~Xt   *    \  exp(X  Trepl)    j  e_Xx  g(t-x)dx 

repl 

Recall   the  Laplace   transformation   is   defined  as    (19) 

£{f(t)}-j     dt   e"St   f(t)    =  f(s).  (2.24) 

0 

The  Laplace  transform  of  Eq.  (2.23)  thus  is 

CO  00  £ 

g(s)   -       dt  e   S     e         +  X  exp(X  T         )    I    dt  e  St         dx  e"  X  g(t-x).  (2.25) 

0  repl  repl 

The   first   term  on  the   right   hand   side   in  the  above   equation   is 


0 

<  T 

< 

repl 
(2.23) 

T 
l 

epl 

< 

t    <  ». 
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FIG.  2.3.  The  numerical  solution  to  the  renewal  equation  for  a 
single  plant  and  different  S  values. 
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j-   "St   -Xt     Tc    -Xt-,  1  ._  „-. 

dt  e    e    =  Lie        }  =  j-^  .  (2.26) 


To  evaluate  the  second  term  in  Eq.  (2.25),  a  variable  is  introduced  such 

that  w=x-T   ,  or  x  =  w  +  T   ,.   Thus, 
repl  repl 


t  t-T 

_»_  !       repl   , 

dx  e   g(t-x)dx  =1    dw  e   exp(-X  T^  J    gCt-W-T^  1).  (2.27) 

repl  0 

A  second  change  in  variables  is  now  introduced,  namely  v  =  t  -  T   ,  so 

repl 

that  the  second  term  on  the  right  hand  side  of  Eq.  (2.25)  becomes 
0  0 


,o°    -s(v+T    )  rv  , 

X   dw  e        '     dw  e  "  g  (v-w) . 


Recall  from  Laplace  transforms  (19) 

Jr  I 

M\      f(u)  g(t-u)duf  =  £{f(u)}I{g(u)}  =  f(u)  g(u),         (2.28) 

where  f(u)  and  g(t-u)  are  some  arbitrary  functions.   The  second  term  on 
the  right  hand  side  of  Eq.  (2.15)  can  now  be  rewritten  as 


X   exp(-s  Tr   )  j   dv  e"SV  j  dw  e~*w  g(v-w)  =  X   exp(-s  T    ) 
0  0 


e  Xw  g(v-w)dw  =  X   exp(-s  T^  ^  ds~Xw\   i(s) 


l^(-sV^i(s)' 


Equation  (2.25)  can  now  be  written  as 
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i(s)  .  _i_  +  x  exp(.s  T^  |M  ,  (2.29) 


which  upon  solving  for  g(s)  gives 


i(s)  =  s  i  X  -  I  exp(-s  T   J  •  <2-3°> 

repl 


The  solution  of  the  renewal  equation  g(t)  can,  in  principle,  be  obtained 
by  taking  the  inverse  Laplace  transform  of  g(s),  i.e., 


g(t)  =  L~l 


Is  +  X  -  Xexp(-s  T   ,) 
repl 


T  >  T   .      (2.31) 
—  repl 


This  inversion  appears  to  be  very  difficult,  or  if  possible,  to  perform 
analytically.   However,  it  is  possible  to  obtain  from  it  the  asymptotic 
solution.   The  asymptotic  solution  for  g(t)  as  t  -*■  ■»  can  now  be  found 
by  application  of  the  final  value  theorem  of  Laplace  transforms  (19) 

lira  f(t)  =  lim  s  1(a)    ,  (2.32) 

t+oo         s+0 

where  f(t)  is  some  arbitrary  function  and  f(s)  =  L{f(t)}. 

Thus, 

lim   .  .    lim   -,  .    lim  ,  S  , ; — - r-  .    (2.33) 

t_  8(t)  ■  8+0  s  g(8)  =  ^0  s  ♦  X  -  Xexp(-sTrepl) 

Application  of  L'Hospital's  rule  to  the  above  limit  yields 

s+0   sv  '        1  +  X  T   ,    1/X  +  T 7 
repl         repl 

Note  that  1/X  is  the  expected  time  the  system  is  operational  and  T 

repl 

is  the  expected  time  the  system  is  down.   This  result  is  reasonable 
because  the  system  should  stabilize  after  a  certain  number  of 
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repair/replacement  time  intervals.   Thus,  the  probability  the  plant  is 

operating  is  equal  to  the  mean  time  to  failure  divided  by  the  total 

cycle  time.  This  asymptotic  solution  for  g(t)  agrees  with  the 

asymptotic  probability  derived  by  Parzen  (17)  for  a  general  distribution 

F   (t)  given  by 
on    °      ' 


„  on    off 


2.3.4  Consideration  of  Multiple  Plant  Spare-Component  Problem 

The  renewal  process  analogy  can  be  extended  to  include  the 
possibility  of  more  than  one  plant  belonging  in  the  spare-component 
pool.   The  N  members  of  the  spare-component  pool  comprise  a  system  of  N 
independent  plants  where  the  probability  the  entire  system  is  operating 
with  original  (or  repaired/replacement)  plant  components  at  time  t, 
denoted  by  g  (t),  is  simply 

gs(t)  =  [g(t)]N,  (2.36) 

where  g(t)  is  the  probabililty  of  a  single  plant  operating  with  an 
original  (or  repaired/replacment)  plant  component  at  time  t.   The 
quantity  g  (t)  is  thus  the  probability  the  spare  component  is  unused  and 
available  should  it  be  needed.   For  the  special  case  of  a  perfect  spare 
(i.e.,  zero  failure  rate),  g(t)  is  given  by  the  renewal  equation 
solution  presented  in  the  previous  section.   Thus,  the  asymptotic 
probability  that  the  entire  system  is  operational  at  time  t  is  given  by 


k    repl7 
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2.3.5  Probability  of  Spare  Availability 

The  probability  that  the  spare  is  available  at  any  time  t  estimated 
by  renewal  theory  is  the  probability  that  all  plants  are  operating  at 
time  t  given  by  Eq.  (2.36).   The  results  of  an  example  problem  are 
graphically  represented  in  Figure  2.4.   A  study  of  the  figure  indicates 
that  the  relative  minimum  value  and  successive  oscillations  in  the  g  (t) 
values  are  sensitive  to  the  product  of  the  component  failure  rate  and 
the  repair/replacement  time  and  insensitive  to  the  number  of  plants  in 
the  spare-component  pool. 

2.3.6  Consideration  of  Multiple  Spares  Problem 

The  renewal  process  analogy  can  also  be  extended  to  include  the 
possibility  a  system  of  N  plants  will  be  operating  when  two  spare 
components  are  available.   The  probability  that  one  spare  will  be 
available  at  the  time  of  a  plant  failure  is  the  sum  of  the  probability 
that  all  plants  in  the  system  are  operating  plus  the  probability  that  at 
most  one  plant  is  down.   The  probability  of  spare  availability  for  the 
case  of  two  spares  is  thus 

gs(t)  =  [g(t)]N  +  N[g(t)]N_1[l-g(t)].  (2.38) 

Substituting  the  analytical  solution  for  g(t)  for  t  <  T     and 

repl 

using  the  asymptotic  solution  to  g(t)  for  large  t,  the  probability  the 
entire  system  is  initially  operating  or  operating  after  the  pool  has 
been  in  place  for  a  long  time  is 


g.(t)  = 


(1  -  <  T 

repl 


At 

(1 

-   e 

*t, 

N- 

1 

1 

t+) 

T         J 
repl' 

(2.39) 


1+A  T   J       ll+X  T   J         1+X  T   J         repl 
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In  general,  for  the  case  of  k  spares  the  probability  that  at  least 
one  spare  will  be  available  is 

8.W  "  I   f?l  tg(t)]N"J  [l-g(t)]3  (2.40) 

j=o  uj 

so  that  the  initial  and  asymptotic  system  operational  probabilities 
become 


Y  [?]  e-^^Cl-e-^V   ,   0  <  t  <  T   . 

jio  UJ  repl 


For  the  remainder  of  this  present  study  only  the  case  of  a  single 
spare-component  pool  is  considered. 

2.3.7  Estimation  of  Fractional  Downtime  Using  Renewal  Theory 

The  fraction  downtime,  denoted  as  £~T,    is  defined  as  the  effective 
fraction  of  the  repair/replacement  time  that  a  plant  is  shut  down  or 
derated  due  to  a  failure.   The  fraction  downtime  can  be  expressed  as 


.   _  effective  plant  shutdown  or  derated  time  due  to  a  failure 

DT  T 

repl 


For  Plan  A,  the  plant  is  shut  down  the  entire  repair/replacement  time 
(i.e.,  fraction  downtime  =  1).   For  Plan  B,  the  plant  may  be  shutdown  or 
derated  only  a  portion  of  the  repair/replacement  time  depending  on  the 
availability  of  the  spare . 

In  this  section  the  fraction  downtime  is  derived  for  the  following 
assumptions : 
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1)  The   removal   time   is    the   same   for   the   failed   component   and   the 
spare.      The   installation  time   is    the   same   for   the   repaired/ 
replacement   component  and   the   spare. 

2)  The  spare's   failure   rate   is   zero    (i.e.,    a  perfect   spare). 

3)  The  spare   is   removed   and   returned   to   storage    (and   thus  made 
available   to   other  plants)   when  the   repaired/replacement 
component    is    installed  and   returned   to  operation. 

Under   the  assumptions   listed  above,    the   fractional  downtime   can  be 
expressed  as 


fraction  of  the  interval 
T  required  to  remove 

the  failed  component,   in- 
stall and  remove  the  spare, 
and  install   the  repair/ 
replacement  component 


Prob- 


spare  ls  not 
available 
at  the  time 
the  failed 
component  is 
removed 


remaining  fraction 
of  the  T 


for  removal  or 
installation 


where 

I  spare  is  not 
available 
at  the  time 
the  failed 
component  is 
removed 


=  1  -  Prob 


spare  is  not 
installed  in 
either  of  the 
other  plants 
in  the  spare 
pool 


spare  has  not 
failed  during 
the  time  it  is  I 
installed  in  a  f 
plant 


The  above  expression  for  spare  availability  can  be  approximated  by 
renewal  theory  as 


Prob 


spare  is  not 
available  at  the 
time  the  failed 
component  is 
removed 


> « i  -  [g(t)]<N-1)(i), 


(2.42) 


where  t  is  the  time  the  failed  component  is  removed,  N  is  the  number  of 
plants  in  the  spare  pool,  and  the  probability  the  spare  has  not  failed 
equals  to  1. 
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Let  T     be  the  maximum  temporary  time  the  spare  can  be  installed 
temp                      J                          r 

and  T   be  the  sum  of  the  failed  component  removal  and  spare 

installation  times.  Thus,  the  repair/replacement  time  interval  is  made 

up  of  a  removal  and  installation  time  at  the  beginning  of  the  interval, 

a  temporary  spare  use  time  in  the  middle,  and  a  removal  and  installation 

time  at  the  end  of  the  interval  (i.e.,  T   ,  =  T     +  2T  ) .   The 

repl    temp     ex 

fraction  downtime  can  now  be  estimated  by 


f   (Njt)  ,  repl  ~  Ttemp   , 

lDTv '  '  T 

repl 


1  -  [g(t+Tex)r'-^^^  .    (2.43) 


For  the  case  of  H  -  1,  f„„  (l,t)  =  (T   ,  -  T    ) /T   .  while  for  the 
DT  repl    temp   repl 

case  of  a  very  large  number  of  plants  (N  *  °°)  ,  f   (»,t)  ♦  1,  i.e.,  the 
spare  is  always  in  use  and  never  available.   For  the  case  of  constant 
g(t)  (i.e.  asymptotic  pool  operation)  the  fractional  downtime  is 


T   i  "  T-  f    i    )N"M  T. 

lmW   .  ^SElt-E  +  ,  -  ^1       La.       (2.44) 

repl  v    repl'      repl 


2.3.8  Sample  Problem  For  a  Five  Component  System 

The  fractional  downtime  as  a  function  of  time  was  calculated  for  a 

system  of  five  plants  with  X  -  .01  yr"  ,  T   ,  =  2  yr.,  and  T 

*      repl     *         temp 

1.833  yr  (T   =2  mo).   The  fractional  downtime  was  calculated  using  the 
numerical  solution  for  g(t)  derived  in  Section  2.3.2  and  the  asymptotic 
solution  for  g(t)  derived  in  Section  2.3.3.   The  results  for  both 
numerical  and  asymptotic  g(t)  are  listed  in  Table  2.1. 
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Table  2.1  Fractional  downtime  estimated  by  renewal 

theory  with  component  failure  rate  = 

.01(yr  ),  repair/replacement  time  = 

1.833  yr,  and  5  plants  in  the  spare 
pool. 


,  fractional  downtime 

time  asymptotic  g(t)        actual  g(t) 


0.0  0.1533  0.0866 

0.1  0.1533  0.0938 

0.2  0.1533  0.1011 

0.4  0.1533  0.1153 

0.6  0.1533  0.1294 

0.8  0.1533  0.1432 

1.0  0.1533  0.1540 

1.2  0.1533  0.1537 

1.4  0.1533  0.1535 

1.6  0.1533  0.1534 

1.8  0.1533  0.1533 

2.0  0.1533  0.1533 

2.2  0.1533  0.1533 

2.4  0.1533  0.1533 

2.6  0.1533  0.1533 

2.8  0.1533  0.1533 

3.0  0.1533  0.1533 

The  time  is  dimensionless  (i.e.,  time  =  t/T   ,) 

repl 

2.4  PC-SPARE  Probability  of  Spare  Availability 

The  results  for  the  fractional  downtime  using  renewal  theory  to 
describe  the  probability  of  spare  availability  can  be  compared  to  the 
method  using  binomial  theory  outlined  in  PC-SPARE  (12).   Recall  that  the 
probability  of  j  failures  in  time  t,  denoted  by  P  (j,t),  is  given  by 

Pf(j,t)  -"\\      e .  (2.45) 

For  the  case  of  constant  repair/replacement  time  and  N  plants  in  the 
spare  pool  the  probability  that  the  spare  is  available  at  time  t  is 
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constant  and  dependent  on  the  number  of  failures  in  the 
repair/replacement  interval  being  analyzed.  The  probability  of  spare 
availability,  denoted  by  P     ,  can  be  expressed  as 


N       (k-1)  failures  in 
I     ProbiT     in  the  other 


k-1 


l(N-l) 


plants 


only  one 
failure  can 
occur  per  plant 


of  any  failed 
■  ProM plant  claiming/ 
the  spare 


Each  plant  is  limited  to  one  failure  because  the  plant  cannot  return  to 

service  (i.e.,  original  component  installed)  by  the  end  of  T   ,  should 
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a  failure  occur.  The  probability  of  a  failure  in  time  t  is 


Pf(l,t)  =  Pf  = 


C*  t) 


1  -At 


1! 


=  X  t  e 


■Xt 


(2.46) 


Thus,  according  to  binomial  theory 


(k-1)  failures  in 

ProMT   ,  in  the  other 
repl  , 

l(N-l)  Plants 


only  one 
failure  can 
occur  per  plant 


t\  ■PnV" 


(N-l)  ! 


(k-l)!(N-k)!    {       repl 


-XT       ,ik-l 


X  T 


"if"'    (l   -   X 


-XT 


N-k 
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repl 


(2.47) 


The  assumption  is  made  that  each  plant  that  fails  has  equal  claim  to  the 
spare,  thus 

of  any  failed     , 


Prob  1  plant  claiming/ 
the  spare 


k  - 


The  probability  of  spare  availability  can  now  be  written  as 


(N-l)  ! 


avail    '  k!(N-k)!  {       repl 


-XT 


k-1 


X  T 
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-XT 


N-k 


1  -  X  T 


repl 
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.  (2.48) 
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The  method  outlined  above  allows  for  a  failure  in  the  interval 
T   ..  and  a  plant  still  has  a  probability  of  obtaining  the  spare.   The 
value  of  P      is  the  average  probability  that  a  plant  would  have  use 
of  the  spare  if  the  plant  were  to  fail  sometime  in  the  interval.   The 
method  says  nothing  about  the  order  of  the  failures  and  hence  the 
decision  as  to  whom  gets  the  spare  is  similar  to  waiting  until  the  end 
of  the  interval  and  all  the  failed  plants  have  "drawn  numbers  from  a 
hat"  to  decide  which  plant  obtains  rights  to  the  spare. 

In  contrast,  renewal  theory  awards  the  spare  to  a  plant  on  the 
"first  come,  first  served"  basis.   Any  failure  in  the  previous  T     of 
the  last  failed  plant  and  the  spare  is  not  available.   For  this  reason 
the  probability  of  availability  estimated  by  PC-SPARE  should  be  larger 
than  the  probability  of  availability  estimated  by  renewal  theory.   For 
the  case  of  a  low  probability  of  a  failure  occurring  in  the  T 
interval  (e.g.,  small  failure  rates,  small  failure  interval,  and  small 
number  of  plants)  the  results  for  probability  of  spare  availability 
estimated  by  PC-SPARE  and  renewal  theory  should  be  approximately  the 
same. 

2.4.1   Sample  Problem  Comparing  Probability  of  Spare 
Availability  Estimates. 

The  following  example  cases  were  calculated  using  both  renewal 

theory  and  PC-SPARE  to  predict  the  probability  of  spare  availability 

(see  Table  2.2).   As  expected,  the  two  methods  for  probability  of  spare 

availability  are  approximately  equal  for  the  cases  of  a  low  probability 

of  a  failure  occurring  in  the  T   ,  interval.   The  PC-SPARE  method 

repl 

estimates  the  average  probability  of  spare  availability  for  all  plants 
in  the  pool.   In  contrast,  the  renewal  theory  method  awards  the  spare  to 
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a  plant  on  the  "first  come,  first  served"  basis  thereby  ignoring 

instances  when  the  spare  is  not  available  for  other  members  of  the  spare 

pool  that  fail  in  the  T   .,  interval. 
repl 

Table  2.2  Probability  of  spare  availability  estimates  using  renewal 
theory  and  PC-SPARE  methods. 


Failure    Repair/replacement 
rate  time 

(yr  )  (yr) 


Number  of       Probability  of  spare 

plants  in  availability 

spare  pool    Renewal  theory    PC-SPARE 


0.01 


0.1 


0.1 


0.5 


0.9 


0.5 


5 

20 

100 

5 

20 

100 

5 

20 

100 

20 


998004 

.996010 

990567 

.981188 

952128 

.905784 

990099 

.980248 

954118 

.909588 

789354 

.610323 

982320 

.964796 

919634 

.843468 

663473 

.411887 

.654579 


.395735 


2.5  Evaluating  Probability  of  Spare  Availability  and  Fractional  Downtime 
Using  Computer  Simulation  Techniques 

The  fractional  downtime  can  be  determined  analytically  only  for  the 

simple  case  of  a  perfect  spare  and  the  condition  that  the  spare  is  only 

used  by  the  "recently  failed"  plant  if  immediately  available.   The 

theoretical  description  that  takes  into  account  all  possible 

combinations  of  pool  life-time,  component  repair  time,  etc.  is  presently 

not  available.   Consequently,  a  computer  code  using  Monte  Carlo 

techniques  was  developed  to  predict  the  fractional  downtime  by 

simulating  random  Poisson  failures. 
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Monte  Carlo  techniques  are  applied  to  dynamic  problems  for  which  a 
closed  form  mathematical  solution  is  difficult  or  impossible  to 
construct  (15).   An  approach  similar  to  that  developed  by  Widawsky  (14) 
and  Emshoff  and  Sisson  (15)  was  considered  where  fixed  and  variable 
time-increment  methods  were  considered.   In  the  fixed  time- increment 
method  (also  known  as  "time  slicing")  time  is  incremented  by  a  small 
constant  interval.   At  the  end  of  each  small  time  step  the  system  is 
checked  to  see  if  any  failures  have  occurred  during  that  interval.  The 
advantage  of  this  method  is  that  a  record  is  not  needed  of  the  sequence 
of  events  because  at  each  time  increment  the  possibility  of  each  event 
is  checked  (14).   The  disadvantages  are  1)  the  time  increments  must  be 
short  in  comparison  to  the  mean  time  between  failures  and  2)  because 
failures  and  repairs  are  randomly  distributed  throughout  time  there  are 
many  intervals  where  no  failures  will  occur. 

In  the  variable  time-increment  method  (also  known  as  "event 
sequencing")  time  is  variably  incremented  to  the  next  predicted  failure. 
The  advantage  of  this  method  is  that  the  time  increments  may  be  long  as 
well  as  short  and  that  a  failure  occurs  only  at  the  beginning  and  end  of 
the  discrete  intervals.   The  disadvantage  is  that  a  record  of  the 
failures  that  occur  must  be  kept  for  the  more  general  conditions  of 
unequal  failure  rates.   Because  of  the  fewer  computations  required  the 
variable  time-increment  method  was  selected  for  use  in  this  computer 
simulation  code. 

The  simulation  is  achieved  by  generating  a  sequence  of  random 
numbers  from  the  appropriate  probability  distributions.  A  Poisson 
process  can  be  simulated  by  recalling  that  the  interarrival  times, 
denoted  by  At,  have  a  cumulative  distribution  function  (CDF)  of  (15) 


F(At)  =  1- 


-XAt 
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(2.49) 


where  A  is  the  individual  component  failure  rate.  For  a  system  of  n 
components  with  the  same  failure  rate  X  operating  at  the  beginning  of 
the  operation  interval  At 

-nAAf- 

"  .  (2.50) 


F(At)  -  1  -  e  ■""*  . 

Thus,  the  time  to  the  next  failure  is 
A(.  „  -Ml-F(At)] 


(2.30 


Because  F(At)  lies  in  the  interval  [0,1],  the  value  for  At.  can  be 

l 

simulated  by  randomly  selecting  a  random  number  u  from  a  uniform 
distribution  over  the  interval  [0,1]  (see  Fig.  2.5). 


F(At.)   u. 


Thus,  the  time  interval  At  between  failures  can  be  simulated  by 
substituting  a  random  number  for  F(At)  in  Eq.  (2.50)  making 


4tt  =[-te(l-Ul)]/nX 


(2.52) 


where  n  is  now  the  total  number  of  components  in  the  system  that  are 
operating  at  the  beginning  of  the  interval  with  a  common  failure  rate  of 
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X.   Successive  operational  intervals  are  analyzed  in  the  same  manner 
until  the  sum  of  the  At  's  are  greater  than  the  maximum  simulation  time 

V 

In  this  manner  a  complete  system  history  over  the  lifetime  of  the 
spare-component  pool  is  obtained.   To  obtain  the  average  or  expected 
characteristics  of  a  particular  pool  scenario,  many  such  simulated 
histories  are  generated  and  the  desired  characteristic  averaged  over  all 
these  histories. 

The  development  of  all  simulation  models  for  this  project  was  based 
on  the  components  failing  with  constant  failure  rates,  thereby  randomly 
distributing  the  failures  according  to  a  Poisson  process.   The  times 
between  failures  for  a  Poisson  process  with  an  infinite  calling  source 
are  distributed  negative  exponentially  characterized  by  a  single 
constant  failure  rate  (17).   In  contrast,  the  systems  considered  in  this 
project  consisted  of  a  small  number  of  components.   The  times  between 
failures  are  distributed  negatively  exponentially  but  each  interval  is 
characterized  by  a  specific  failure  rate  depending  on  the  number  of 
components  operational  at  the  beginning  of  that  interval.   Thus,  the 
time  intervals  between  failures  have  different  negative  exponential 
distributions. 

The  cumulative  distribution  function  F(t)  is  defined  as  the 
probability  a  failure  will  occur  in  the  interval  (0,t).  Thus, 

a  single  failure  in    r 
F(t)  -  Probj  =    dt'  f(t'),       (2.53) 

the  interval  (0,t)     > 

where  f(t')dt'  is  the  associated  distribution  density  function.   This 
density  function  is  the  probability  of  a  component  failing  in  dt '  about 
t\ 
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An  important  property  of  the  Poisson  distribution  is  the 
"forgetfulness  property".  Consider  a  plant  component  that  is  Poisson 
with  mean  rate  A  of  two  failures  per  day.   If  the  component  has 
experienced  no  failures  for  the  previous  t  hours  of  operation,  what  is 
the  probability  it  will  fail  in  the  next  t  hours,  i.e.  what  is 


ProbfT  >  t   +  t  I  T  >  t  }  (2.54) 

P  P 

where  T  is  the  time  from  t  to  the  time  of  the  next  failure.  Recall 
P 

that  the  conditional  probability  of  A  given  that  B  has  occurred  is  (20) 
Prob(A|B)  =  P(Ap^  B>  .  (2.55) 


Thus,  Eq.  (2.54)  can  be  written  in  the  form 

Prob(T  >  t  +  t  and  T  >  t  } 
Probd  >  tp  ♦  t  |  t  >  tp> progtT  ;  ■  , E_  .       (2.56) 


A  closer  analysis  of  Prob(T  >  t   +  t  and  T  >  t  )  shows  that 

P  P 

Prob(T  >  t  +  t}  necessarily  includes  Prob{T  >  t  },  therefore,  Eq. 
P  P 

(2.56)  becomes 


Prob{T  >  t   +  t} 
Prob{T  >  tp  ♦  t  |  I  >  tp>  =   Prob(T  >Pt  ;   •         (2.57) 


Recalling  that  the  intervals  between  failures  are  negative  exponentially 
distributed,  yields  for  Eq.  (2.57) 


-X(tp+t) 


At 


ProbfT  >  t   |  T  >  t  }  =  - — r- =  e"  Z    .  (2.58) 

P        P      -At 


The  results  of  Eq.  (2.58)  show  that  waiting  time  t  has  no  effect  on  the 

P 

probability  of  a  failure  in  the  next  t  interval.   This  independence 
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between  two  intervals  of  time  is  the  so-called  "forgetfulness  property" 
Thus,  for  each  operation  interval  simulated,  the  only  concern  need  be 
the  number  of  components  operational  at  the  beginning  of  that  interval 
and  no  consideration  need  be  given  the  previous  failure  history.  Thus, 

D,-„k/of  a  failure  I  a  failure  occurred  \       _  .  /of  a  failure! 
Prob\in  (Tl,T1+t)  |  at  T±  /  =  Prob\ln  (O.t)    /• 

2.5.1  Case  of  a  System  Without  a  Spare 

Under  Plan  A,  one  must  estimate  the  probability  of  a  failure 
following  the  first  failure  in  a  system  where  the  repair/replacement 
time  for  the  failed  component  is  a  fixed  time  T     (see  Fig.  2.6). 

t' 


T  =0  T  J  L       T  +T   ,        tlme 

°  dt'  rep 

FIG.  2.6.  A  graphical  representation  of  a  failure  in  dt'  about  t' 
following  the  first  failure  at  T  in  a  system  with  a  fixed  repair 
time. 

For  this  discussion  assume  at  time  T  there  are  N  identical 

o 

components  operating  with  identical  failure  rates  A.  The  system 
characteristic  failure  rate,  denoted  by  rate  X,  is  defined  as  the  system 
failure  rate  in  a  very  small  period  of  time  dt'  and  is  the  sum  of  all 
the  individual  component  failure  rates  operating  at  the  beginning  of  the 
interval  dt'.  Thus,  X   for  the  time  intervals  following  T  is 

0  <  t'  <  T 

"  "Pl  (2.59) 

,   T     <  t'  <  »  , 
repl      — 
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where  n  is  the  number  of  plants  operating  at  the  beginning  of  the 
interval  dt' .  Thus, 


Prob{of  a  failure  in  dt '  about  t ' }  =  A  dt '  =  NA  dt ' 


(2.60) 


The  density  function  for  the  time  between  failures  can  be  written  as 
f(t')dt'  =  Probj 


, /of  a  failure  in  dt1  j  a  failure 


\about  t' 


occurred  at  T 


} 


N 


[  Prob 


i=l 


I  of  a  failure  of  the  i-th 
j component  in  dt*  about  t' 


the  i-th  component 
survives  until 
at  least  t' 


(2.61) 


The  probability  a  component  survives  to  time  t'  (measured  from  T.)  is 
e    .   Thus,  the  above  equation  can  be  rewritten  as 


f(t')dt'  =  {NX  dt'He~Xlt'e"X2t'.  .  .  e  N   }  (2.62) 

where  N  is  the  number  of  operating  elements  at  t'.   The  density  function 
for  the  time  between  failures  is  defined  for  two  values  of  A  thereby 
resulting  in  the  density  function  defined  for  two  intervals  of  interest. 


Case   1:      0  <  t'    <  T 


repl 


f(t')dt'   =   {(N-l)A  dt'He 


■X1t'      -A2t' 


-Vic\ 


(N-l)A   e-(N-lH   f    dfc, 


(2.63) 


Case   2:      T        ,    <   t'    < 
repl  — 


The  N-th  plant  becomes  operational  only  at  time  T     after  T.  and, 

repl  1 


therefore,    only  has   to  survive    (t'    -  T        ,)    making 

repl 


f(t')dt'    =   {NA  dt'He 


■Xyt'  -A^' 


-A        t'      -A.,(t'-T        ,) 
N-l  N  repl    , 
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-  NX  e-(N"1)X  C'  ."X(t''T«P1) 


XT 
.,,   -NX  t'    repl  ,.  ... 

NX  e       e    r  .  (2.64) 


Thus,  the  density  function  for  the  time  between  the  first  and  second 
failure  is  defined  for  the  two  ranges  of  t'  as 

[(N-nx  e-<N-i)*f  dt-       ,     o  <  t<  <T 


««'>*'    ■  1  _■«,.      X  T™!  r  rePl  (2.65) 

repl 


IBX  e-NXt'    e        repl  df      , 


If  a  substitution  of  variables  is  made  with  B  =  XT   , .  T  =  Xt,  and 

repl 

dxf  =  Xdt',  the  above  can  be  written  as 

f(N-l)  g-^-1^'  dTi   o  <  t  <  B 
t(l')   df  =  (2.66) 

N  e  e     dx'  B  <  X   <   » 


The  cumulative  distribution  function  for  this  density  function  is 

4t  f«  -(N-l)Ax 

C  ■  l-ev  ,0<At<6 

F(4t)    =  dt'    f(T')    =  (2.67) 

J0  [(1  -  eB  e"NAT)        ,      B   <  At   <  -   . 

The  results  for  F(At)  satisfies  the  usual  condition  F(0)=0,  F(=°)  =  1  and 
is  continuous  at  At  =  B. 


2.5.2  General  Description  of  Simulation  Models 

A  similar  procedure  for  finding  f(T')  dT'  and  F(t)  is  followed  for 
the  situation  in  which  a  spare  is  added  to  the  system  (Plan  B) .   The 
density  function  to  be  derived  is  a  function  of  the  number  of 
operational  components  at  the  beginning  of  each  time  interval  of 
interest.   The  time  interval  of  interest  is  the  time  interval  between 
the  last  simulated  failure  and  the  next  failure  to  be  simulated. 
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The  terminology  of  the  "k-th"  plant  was  selected  to  mean  the  plant 
for  which  the  probability  of  spare  availability  or  fractional  downtime 
is  to  be  estimated.   In  addition,  the  repair/replacement  time,  denoted 
by  3,  is  a  standard  dimensionless  time  unit.   The  interaction  between 
the  other  (N-l)  plants  and  k-th  plant  is  simulated  to  obtain  an  estimate 
of  the  probability  of  spare  availability  and  fractional  downtime  for  the 
k-th  plant.   To  determine  this  interaction  three  models  were  developed 
to  simulate  plant  failure  histories. 

Model  1  is  the  simplest  simulation  performed.   The  major 
assumptions  are: 

1)  The  failed  component,  spare,  and  repaired/replacement 
component  are  removed  or  installed  instantaneously. 

2)  The  spare's  failure  rate  is  zero  (i.e.  a  perfect  spare). 

3)  The  spare  is  immediately  installed  in  a  failed  plant  when  it 
becomes  available. 

Because  the  spare  has  a  zero  failure  rate,  the  system  appears  to  have 
only  as  many  components  as  there  are  plants.   Also,  the  spare  is 
available  if  all  plants  are  operating  at  the  time  of  the  k-th  plant 
failure.   The  assumptions  used  in  this  model  are  the  same  assumptions 
used  in  the  renewal  theory  approximation  to  the  spare-component  pool 
problem.   As  a  result,  the  probability  of  spare  availability  estimated 
by  this  model  should  agree  with  the  probability  of  spare  availability 
estimated  by  renewal  theory.   However,  the  fractional  downtime  estimates 
will  not  agree  because  renewal  theory  does  not  consider  the  possibility 
that  the  spare  may  become  available  prior  to  the  end  of  the 
repair/replacement  interval. 
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In  Model  2  the  major  assumptions  are: 

1)  The  failed  component,  spare,  and  repair /replacement  component 
are  removed  and  installed  instantaneously. 

2)  The  spare's  failure  rate  is  zero  (i.e.,  perfect  spare). 

3)  If  the  spare  is  not  immediately  available  at  the  time  of 
failure  it  is  never  used  to  cover  any  portion  of  the  downtime 
for  that  failure. 

Once  again,  because  the  spare  has  no  failure  rate,  the  system  appears  to 

have  only  as  many  components  as  there  are  plants. 

The  availability  of  the  spare  to  the  k-th  plant  upon  failure  is 

determined  by  the  last  time  the  spare  was  used.  The  spare  is  available 

if  there  were  no  failures  within  6  (where  6  =  XT   , )  of  the  k-th 

repl 

plant's  failure,   (see  Fig.  2.7). 


Tk        '   Tk  6  Tk 

FIG.  2.7.  Representation  of  a  failure  sequence  where  t  is  time 
of  last  failure  before  the  k-th  plant  failure. 

If  the  spare  is  not  used  by  any  of  the  plants  that  fail  in  the 
interval,  the  spare  is  available.  The  spare  is  also  available  if  the 

last  time  the  spare  is  used  is  outside  of  the  Interval  (t,  -6,  t,  )  and 

k      k 

the  difference  between  the  failures  that  occur  in  the  previous  interval 
are  within  8  of  the  last  time  the  spare  was  used  (see  Fig.  2.8). 


used      k         L     used 


FIG.  2.8.  Representation  of  a  failure  sequence  where  T  is  the 
time  of  the  last  failure  and  t  is  the  last  time 
the  spare  was  used  prior  to  thesfc-th  plant  failure. 
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Unlike  Model  1,  the  spare  may  be  available  even  when  failures  occur 

in  the  interval  (t,  -  B,  T.) .  However,  for  the  spare  to  be  available 

t    ,  cannot  fail  in  the  interval  (t.  -  S.  t,  ) .  The  results  from  this 
used  k      k 

model  should  show  a  larger  probability  of  spare  availability  than  Model 
1  or  renewal  theory  because  the  spare  may  be  in  storage  at  t  with  all 
except  one  plant  failed. 

For  Model  3  the  major  assumptions  are: 

1)  The  failed  component,  spare,  and  repair/replacement  component 
are  removed  and  installed  instantaneously. 

2)  The  spare's  failure  rate  is  equal  to  the  failure  rate  of  the 
other  similar  components  when  installed  in  a  plant  system. 

3)  The  spare  is  immediately  installed  in  a  failed  plant  when  it 
becomes  available. 

Model  3  is  identical  to  Model  1  except  the  spare  component  can  fail 
and  the  fraction  downtime  can  now  assume  any  value  over  the  interval 
[0,1];  that  is, the  spare  contributes  to  the  total  number  of  components 
that  can  fail  when  installed  in  the  system.   As  a  result,  failures  can 
now  occur  in  the  interval  (t  -  6,  x  )  and  the  spare  still  be  available 
for  use  by  the  k-th  plant. 

A  system  with  a  finite  spare  failure  probability  should  have  a 
larger  fractional  downtime.   In  contrast,  a  system  that  uses  the  spare 
when  available  should  have  a  smaller  fractional  downtime.   Unlike  the 
previous  models,  Model  3  is  limited  to  evaluating  the  fractional 
downtime  and  not  the  probability  of  spare  availablility. 

2.5.3  Perfect  Spare  Used  When  Available,  Model  1 

In  this  section  a  mathematical  model  is  derived  to  simulate  the 
time  between  component  failures  for  simulation  Model  1.   In  this  model 
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the  perfect  spare  case  is  considered  with  instantaneous  removal  and 
installation  times.  The  spare  is  installed  in  a  failed  plant  when  it 
becomes  available,  even  if  for  a  very  short  period  of  time. 

Two  Plant  Problem 

Consider  first  a  system  of  two  plants  and  one  spare  component.   The 
system  characteristic  failure  rate,  denoted  by  X,  is  defined  as  the  sum 
of  the  individual  component  failure  rates  of  the  plants  operational  at 
the  beginning  of  the  time  interval.   The  number  of  plants  in  the  pool, 
denoted  by  N  (here  equal  to  2) ,  will  remain  constant  throughout  the  pool 
lifetime. 

The  time  between  the  first  failure  and  the  simulation  starting  time 
is  a  waiting  time.   If  all  plants  are  operational  at  time  0,  the  number 
of  plants  operating  at  the  beginning  of  the  interval  is  N.   The  density 
function  for  the  first  interval  can  be  expressed  by 


f(T*)dt'  =  N  dt'(e   L  .  .  .  e  N)  =  N  dT'  e"NT ' 


(2.68) 


The  cumulative  distribution  function  F(At; 

1  is 

f(at; 

At           At 

I  =  J  dT'  f(T')  =  J   N 

0            0 

JT,   -NT' 
dT   e 

=  1 

-1 
-  e 

Mt 

(2, 

.69) 

Thus,  the 

random  number  u  is  substituted 

for 

f(At; 

)  and 

the  first 

simulated 

time  interval  is 

At   =  -£n(l- 

■UjJ/N    , 

0 

<  At 

<   CO 

I  - 

(2. 

70) 

The  s 

;econd  failure  time  is 

dependent 

upon 

the 

first 

failure. 

I 

'or  a 

period  of 

time,  S,  following  the  first  faJ 

.lure 

the 

number  of  operatd 

.onal 

components 

1  will  be  decreased  by 

one,  or 

0  <  T  <^ 

8  <  T  < 
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(2.71) 


and  the  density  function  is 

((N-l)dT'  e^N-D1'     t        0  <  T  <  g 

f<T'>  dT'  "         _NT.   a  (2-72> 

N  dx1  e     e        ,   S  <  T  <  » 


Thus, 

fl  -  e-(N"1)AT  ,   0  <  Ax  <  | 

F(AT)  "  ,     8  -NAT        ,   ,  <2-") 

1-ee  ,   8  <  At  <  »  . 


The  value  of  u-  can  now  be  compared  to  the  probability  of  a  single 
failure  in  the  interval  (t1,  t'  +  8)  to  determine  which  form  of  F(At)  to 
invert  to  obtain  the  simulated  time  interval.  The  result  is  two 
potential  simulation  equations.   The  simulation  equation  is  determined 
by  comparing  u.  to  the  CDF  evaluated  at  increments  of  8.  There  are  two 
cases  to  be  considered. 

Case  1:  u„  <  F(8) 

The  second  time  interval  for  this  case  is  given  by 

At2  =  -[Jtn(l-u2)]/(N-l)  (2.74) 

and 

T2  =  Tl  +  AT2  =  Tl  "  [*n(l-u2)]/(N-l)  .  (2.75) 

Case  2:   F(B)  <  u   <  1.0 

The  second  time  interval  for  this  case  is  given  by 

At2  =  [6  -  *n(l-u2)]/N  (2.76) 

and 
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T2  "  T!  +  4t2  -  Tj  +  [B  -  An  (l-u2)]/N  .  (2.77) 

A  summary  of  the  results  for  the  second  time  Interval  is, 

T  -[An(l-u  )]/(N-l)     ,   0  <  u  <  F(B) 

T2  "  I  (2.78) 

ITj  +  [B  -  An(l-u2)]/N   ,   F(B)  <   u2  <  1.0  . 

The  third  failure  is  also  dependent  upon  the  first  two  failures 
(see  Fig.  2.9).  The  value  of  6     is  defined  as 

«x  »  t   -  (t4  -  t^j)  =  B  -  AT1_1  (2.79) 

where  i  ■  2,3,...(N-1)  and  8,  =  6. 


^r 


T=0     T  T     +|  |+      T+B        T.+B  tlme 

Z     dt'        1  2 

FIG.  2.9.  A  time  line  illustrating  a  possible  configuration 
to  be  considered. 


The  limit  on  i  is  (N-l)  because  at  least  one  plant  must  be 
operational  for  any  failure  to  occur.   If  there  are  no  operational 
plants  then  a  period  of  time  must  elapse  until  at  least  one  plant 
returns  to  operation.  Thus, 

Sj  -  0  -  (tj  -  Tx)  -  ft  -  At2  .  (2.80) 

For  the  case  of  «2  >  0,  a  failure  cannot  occur  in  the  system  until  at 
least  one  component  has  returned  to  service  at  time  I,  +  B  (i.e.,  a 
repair/  replacement  component  is  installed).   The  waiting  time,  denoted 
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by  I  ,  is  given  by  <S„.  The  next  failure  is  then  simulated.  Starting 
from  time  (t   +  6),  follow  the  same  technique  as  outlined  for  the  second 
failure,  and  compare  u,  to  F(8). 

For  the  case  that  S     <   0,  the  plant  that  failed  at  T  has  returned 
to  service  before  t  occurred.  The  density  function  for  the  third 
interval  is  the  same  form  as  given  in  Eq.  (2.72)  where  now 


f(t')  dt'  = 


[(N-l)dT'  e"(N-1)T' 
N  dx-  e-NT' 


0  <  T  <   & 


'  1 


(2.81) 


and 


\l   .  „-(H-1)At 

F(AT)=   ,     \     -NAT 

1  -  e   e 


0  <  At  <  & 


6   <  At  < 


(2.82) 


The  value  of  u  is  then  compared  to  F(t)  to  obtain  the  third  time 
interval 


At, 


Thus, 


|-  In   (l-u3)/(N-l) 
f[«j  -  te(l-u3)]/N 


T2  -  te(l-u3)/(N-l) 
lt2  ♦  [«j  -  tn(l-u3)]/N 


F((51)  <«3<  1.0 


0  <  u3  <  F(«1) 
F(<51)  <Uj<  1.0 


(2.83) 


(2.84) 


General  N-plant  Problem 

The  time  of  the  i-th  failure  can  now  be  determined.   For  the 
general  case  of  N  plants  in  the  spare  pool  and  P  previous  failures  in  a 
time  interval  equal  to  8  just  prior  to  the  failure  of  the  k-th  plant 
component,  the  density  function  for  the  j-th  time  interval  is 


46 


f.(t')  dr'  =  a  exp(-a  t')  exp  \    V  S_  . 
J  i   P+3-m 

m— 1 


■  «P+3_.,  «  f  <  5p.j+2  , 


(2.85) 


where  6 


t  -  6  -  J  At±,  a  =  N-P-2+j,  6p+2  -  0,  $l   -  S,  and  6Q 


Thus,  the  general  form  of  the  CDF  is 

i 

F 


.(A^)  =  1  -  exp(-a  4Tj)  exp^  «p+3_mj ,  6p+3_.  <  At.  <  6p_.+r  (2.86) 


Because  of  the  finite  number  of  plants  in  the  pool  the  maximum 
number  of  prior  failures  in  the  interval  6  is  (N-2) .   If  there  are  (N-l) 
failures,  a  waiting  time  must  be  evaluated  by  subtracting  6_  .  from  the 
6  values  in  Eqs.  (2.85)  and  (2.86).   If  there  are  no  prior  failures 
(e.g.,  P-0),  the  CDF  is 


-NAT. 


F  (At±)  -  1  -  e 


(2.87) 


The  At  for  the  i-th  interval  is  now  found  by 

P  =  0 
j 

K 


-ln(l-u  )/N 


*t-i 


'  I      (5Dj-   -  in[l-u.]l 

",   P+3-m         l 
m-1  J 


[N-P-2+K] 


P  >  0 


such  that  PjWp.W   i  uj    <  FjC«p_k+2)    where   F«P+2>    "   I- 


(2.88) 


Because  the  spare  is  installed  when  it  becomes  available,  the 
downtime  for  the  k-th  plant  depends  only  upon  the  time  of  the  last 
failure.   Consider  the  case  where  t   is  the  time  at  which  the  evaluation 
is  made.   If  the  time  of  the  previous  failure  is  greater  than  B  prior  to 
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t  then  the  downtime  is  0.   If  the  last  failure  has  occurred  within  8  of 
Tk>  say  t._,  then  the  downtime  is  8-(x,  -  T,,).  Thus,  the  downtime  SI 
is  given  by 


(2.89) 

6  -  (tr  -  tlf)   ,   tK  -  tlf  <  6  . 

At  the  end  of  M  simulated  pool  histories  the  fractional  downtime, 
denoted  by  fDT>  is  determined  by  comparing  the  total  simulated  downtime 
to  the  maximum  downtime  MS,  or 

M 


1     S1J 


fDT  "   MB—  •  <2-90> 


The  probability  of  spare  availability  estimated  by  this  model  can 
be  compared  to  the  renewal  theory  results.   The  spare  is  available  at  x, 
if  the  last  failure  occurred  outside  the  interval  (x,  -  8,  x,  ) .   A  value 
of  1  is  assigned  to  a  counter  when  the  last  simulated  failure  prior  to 
Tfe  occurs  in  the  interval  (x,  -  B,  x,  ) .  A  value  of  0  is  assigned  if  the 
last  simulated  failure  occurs  outside  the  interval  (x,  -  B,  x,  )  .   The 
counters  are  then  summed  and  divided  by  the  number  of  simulation.   This 
yields  an  estimate  for  the  probability  of  spare  availability  at  time  x,  , 
denoted  by  P    iiCxt)  and  represented  as 

1  M 

PavailCTk>  ,1"ll  («™t«)j  •  (2-91) 

j=l 

2.5.4  Perfect  Spare  Which  is  Used  Only  if  Immediately  Available, 
Model  2 

In  this  section  a  mathematical  model  is  derived  to  simulate  the 

time  between  component  failures  for  Model  2.   The  same  assumptions  as 
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for  Model  1  are  used,  except  that  if  the  spare  is  not  immediately 
available  at  the  time  of  failure  it  is  never  used  to  cover  any  portion 
of  the  downtime  for  that  failure. 

Because  the  spare  has  no  failure  rate,  the  spare  does  not 
contribute  to  the  system  failure  rate.  Thus,  the  time  intervals  and 
failure  times  can  be  simulated  using  the  method  outlined  for  Model  1. 
In  addition  to  the  failure  times,  the  times  the  spare  is  used  must  be 
estimated.   The  first  failure  simulated  would  use  the  spare  during  the 
fixed  downtime  6.   The  remaining  failures  up  to  the  maximum  evaluation 
time  (pool  lifetime)  would  use  the  spare  only  if  it  were  available  at 
the  time  of  failure.   Thus,  at  each  simulated  failure,  the  previous  time 
interval  of  length  S  is  analyzed  to  determine  if  any  other  failures  have 
occurred  in  the  time  interval.   If  no  other  failures  have  occurred  in 
the  time  interval,  the  spare  is  used  at  the  newly  simulated  failure. 

In  contrast  to  Model  1,  the  spare  can  be  available  even  if  failures 

occur  in  the  interval  (x  -  S,  t  ).   A  counter,  denoted  by  S2,  is  used 

to  determine  the  number  of  times  t    .  occurs  in  the  interval  (t,  -  S, 

used  k 

Tk) .   If  Tused  occurs  in  the  interval  (x  -  g,  x  )  the  spare  is  not 
available  at  x  .   The  probability  that  the  spare  is  available  at  x   is 
given  by 

M 

Pavail(V    'l'*£    S2j    '  <2'92> 

J-l 

where  M  is  the  total  number  of  simulations  run.   Because  the  spare  is 

either  used  if  immediately  available  or  never  used,  the  spare  is  either 

used  the  entire  repair/replacement  time  or  not  at  all.   Thus,  the 

fractional  downtime  is 
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MB 

f. 


M 

J 


M 


DT W   "I-sJ,8V  <2-93> 


Thus,  the  Model  2  probability  of  spare  availability  and  fractional 
downtime  are  equal. 

2.5.5   Imperfect  Spare  Which  is  Qsed  when  Available,  Model  3 

In  this  section  a  mathematical  model  is  derived  to  simulate  the 
time  between  component  failures  for  simulation  Model  3.   Unlike  Model  1, 
the  spare  failure  rate  is  now  equal  to  the  failure  rate  of  other  similar 
components  when  installed  in  a  plant. 

As  in  Models  1  and  2,  the  system  characteristic  failure  rate  is  the 
sum  of  the  individual  component  failure  rates  for  the  plants  operational 
at  the  beginning  of  the  time  interval.   However,  in  this  case  the  spare 
component  contributes  to  the  system's  failure  rate  when  it  is  installed 
in  a  plant. 

Two-plant  Problem 

Consider  a  system  of  two  plants  and  a  spare  component.   The  time 
between  the  first  failure  and  the  simulation  starting  time  is  a  "waiting 
time".   If  all  plants  are  operating  at  time  0,  the  density  function  for 
the  first  time  interval  can  be  expressed  as 


f(T')  dr'  =  N  e"NT'dT'  .  (2.94) 


Thus, 

At  At 


-NAT 

-pIT  "  ,   t 

e    dT'  =  1  - 
0  "0 


f  f       Mt'  "NAT 

F(At)  =  I    dT1  f(T')  =   N  e"WT  dT'  =  1  -  e    .  (2.95) 


50 

The  probability  of  a  single  failure  in  the  interval  (0,»)  is  1.  The 
interval  simulation  equation  is  found  by  substituting  u.  for  F(At)  and 
solving  for  At  in  Eq.  (2.95) 

tj  -  ATj  =  -  dn(l-Ul)/N,     0  <  At j  <  °°  .        (2.96) 

Because  the  spare  now  has  a  failure  rate,  the  number  of  components 

operating  at  the  beginning  of  the  second  failure  interval  is  still  N 

(see  Fig.  2.10).  The  density  function  for  the  second  time  interval  is 

identical  to  that  of  the  first  time  interval,  making 

-NAT 
F(At2)  =■  1  -  e         0  <  At2  <   »  (2.97) 

and 

At2  -  -  ln(l-u2)/N  .  (2.98) 

The  value  of  T.  is 

T2  *  Tl  +  AT2  =  Tl  -[*n<l-u2)]/N  •  <2'") 


T=0 


time 


dT' 


FIG.  2.10.   A  time  line  illustrating  a  possible  configuration  for 
the  second  failure. 


The  density  function  for  the  third  time  interval  is  dependent  upon 
the  value  of  At2.   If  At2  is  greater  than  8,  the  spare  is  available  at 
t2>  otherwise  the  spare  is  not  available  at  the  start  of  At  .   Thus, 
there  are  two  possible  cases  to  be  determined  by  the  value  of  <5  ,  where 
(see  Fig.  2.11). 

<Sl   =  Tl   +  B~T2  =   e~  AT2   ■  (2.100) 


Case  1:   6  <  0 
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For  this  case  the  density  function  for  the  third  time  interval  and 
CDF  are  the  same  as  for  the  first  and  second  failure,  that  is, 

-NAt„ 


F(&t3)  -  1  -  e 


0  <  At  < 


(2.101) 


Case  2:   5  >  0 


Unlike  the  previous  models,  the  spare  now  has  a  finite  failure 
rate.  Therefore,  one  component  is  operating  at  the  beginning  of  the 
interval  and  the  density  function  for  the  third  time  interval  is  now 
defined  over  two  intervals  of  t'.  Thus,  the  density  function  for  the 
third  time  interval  is  given  by 


f(T')  dT"  = 


»-»  e"(N-1)T'  dT- 

I  H  e"NT'  I1     df 


,   0  <  At3  <  i 
,      &l   <  At  <  » 


(2.102) 


Thus, 


F(At3)  =  ■ 


-(N-I)At, 


I   -NAT3 
1  -  e   e 


0  <  At3  <  6 


I,  <  AT3  < 


(2.103) 


At, 


1_i 


T=0 


dT' 


T1+e 


FIG.  2.11.   A  time  line  illustrating  a  possible  configuration 
for  the  third  failure. 


The  value  of  u3  is  then  compared  to  F(6  )  to  find  the  appropriate 
expression  for  F(At3).   The  interval  simulation  equation  is  found  by 
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substituting  u,  for  F(At  )  and  solving  for  At,.  Thus, 


AT3  = 


-  in(l-u3)/(N-l) 
[il  -   in(l-u3)]/N 


0  <  u3  <  F(6L) 
FCiSj)  <  u3  <  1.0  . 


(2.104) 


The  value  of  T,   is  determined  by  S,  which  can  be  expressed  for  two 
possible  cases. 


Case  1:  6     <  0 


T3  "  T2  "  **<l-»3)Al 


(2.105) 


Case  2:  6     >  0 


T3" 


t2  -  in(l-u3)/(N-l)     ,  0  <  u3  <  F(6j) 
t2  +  [«j  -  Jtn(l-u3)]/N  ,  FfSp  <  u3  <  1.0 


(2.106) 


The  density  function  for  At  is  dependent  upon  At,  and  At   (see 
Fig.  2.12).   If  At2  is  greater  than  6  and  At3  is  less  than  8,  the 
interval  is  analyzed  in  the  same  manner  as  described  for  At  .   If  At,  is 
greater  than  6,  the  interval  is  determined  as  described  for  At  or  At  . 
A  different  case  must  be  analyzed  when  At   +  At   <  B. 


I-i-At^I*-  At2 


t=0 


fr 


r    +|  |*  t  +6  t„+B 

J      dT'    *      l 


FIG.  2.12.   A  time  line  illustrating  a  possible  configuration  for 
the  fourth  failure. 


The  6  values  for  the  fourth  failure  are 
S, 


1 


-  At, 


62  =  tx  +  8  -  t3  =  6  -  At3  -  At2  . 


(2.107) 
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The  density  function  for  the  fourth  time  interval  is  now  defined  over 
three  intervals  of  t'.   There  are  three  possible  combinations. 


1.   8,  <  0  and  6.  <  0 


For  this  case 


NAt. 


F(At4)  =  1  -  e 


and 


0  <  At,  < 
4  — 


(2.108) 


At4  =  -  j>n(l-u  )/N 


(2.109) 


2.  6  <  0  and  «,  >  0 


This  case  is  similar  to  the  At,  interval  where  the  density  function 
for  the  fourth  time  interval  is 


f(T')   dT' 


J(N-l)  e~(N  1)T'dT'  ,   0  <  At4  <  Sl 

]„  -Nt'     V 
N  e     e   dT'     ,   6   <  At,  <  » 


(2.110) 


and 


F(At4) 


-(N-I)At. 


1  -  e 


6,   -NAT, 

i     1      * 

1  -  e   e 


0  <  At,  <  6, 
4  —  1 


>1  <  AT4  ± 


(2.111) 


The  value  of  u4  is  then  compared  to  F(6  ).  The  interval  simulatio 
equation  is  found  by  substituting  u,  for  F(AO  and  solving  for  At,, 


thus 


AT4  = 


-  Zn(l-u4)/(N-l) 
[Sj  -  Zn(l-u4)]/N 


0  <  u4  <  F(6:) 
F(61  <  u4  <  1.0 


(2.112) 
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3.  <52  >  0  and  6  >  0 

For  this  case  all  system  components  and  the  spare  have  failed.   As 
in  previous  models,  a  waiting  time  is  incurred  until  a  plant  component 
is  operational.   The  waiting  time  is  given  by 

Tw  =  Tl  +  B  _  T3  •  (2-113) 

Therefore,  the  next  simulation  is  started  at  time  x„  =  T.  +  t  ■  x,  +  S. 

3    3    w    1 

The  next  simulated  time  interval  uses  the  procedure  outlined  above  when 
S,  <  0  and  6      >   0. 

The  density  function  for  Ax.  is  dependent  upon  the  previous  Ax's 

for  a  period  of  I   .  to  x .   ,  where  N  is  the  number  of  plants  in  the 

g       1-1       1-N 

pool.   If  Z      6   <  6,  then  all  components  in  the  system  including  the 
1*1 

spare  have  failed  by  x    .   The  system  must  then  wait  until  T.  .  +  6 
1*"1  l-l    n 

for  at  least  one  component  to. return  to  operation  before  the  next 
failure  time  can  be  simulated.   The  sum  of  the  individual  <5's  must, 
therefore,  be  determined  at  each  x  prior  to  x.   with  the  smallest  x 
being  x    .   The  x  value  at  which  the  sum  exceeds  g  is  the  starting 
index  for  determining  the  number  of  components  operating. 

The  time  of  the  i-th  failure  can  now  be  determined.   For  the 
general  case  of  N  plants  in  the  spare  pool  and  for  P  previous  failures 
(P  >  1)  in  an  interval  of  time  equal  to  6  just  prior  to  the  failure  of 
the  k-th  plant  component,  the  density  function  for  the  j-th  time 
interval  is 


J 


fjCx'Jdx'    =  a   exp    (-ax')    expj    [      &?+2_\ 


«P+2-j   ±T'   iWj    '  (2"U4> 
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where  8,  ■ 


J  At.,  a  =  (N-P+j-1),  <S    =  0,  and  6 


i-1 


P+l 


Thus,  the  general  form  of  the  CDF  is 

f.        1 
F.(At.)  =  1  -  exp(-ax)  «JJ  4     |,  6p_l+2  <  &T£  <  6p_£+1  .  (2.115) 

Because  of  the  finite  number  of  plants  in  the  pool  the  maximum 
number  of  prior  failures  in  the  interval  S  is  (N-l) .   In  addition, 
because  the  spare  component  contributes  to  the  system  characteristic 
failure  rate,  the  general  form  for  the  density  function  and  CDF  apply 
only  to  interval  equations  where  at  least  two  prior  failures  have 
occurred  in  the  interval  (T._,  -  S,,  t   ).   If  there  are  no  prior 
failures  (e.g.,  P  =  0),  the  CDF  is 
-NAT, 


Fj(Ari)  =  1 


(2.116) 


The  value  of  At.  is  found  by  substituting  u.  for  F.(At  )  in 
i  j      r     I' 

either  Eq.  (2.116)  or  (2.115)  and  solving  for  At  giving 


*t,- 


- 

in(l 

-Uj)/N 

K 

Z 

m-1 

6P- 

2+m  "   in(1- 

V 

N   - 

P  -    1   +  K 

P  =  0,1 


y6p-K+2'  ^jiy^-K.! 


(2.117) 
)  and  P>1. 


Because  the  spare  is  installed  when  it  becomes  available,  the 
downtime  for  the  k-th  plant  depends  only  upon  the  time  of  the  last 
failure.   Therefore,  the  fractional  downtime  and  probability  of  spare 
availability  are  calculated  in  the  same  manner  as  used  by  Model  1. 
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2.5.6  Error  Estimates  for  Probability  of  Spare  Availability 

In  this  section  a  method  is  derived  to  evaluate  the  error 
associated  with  the  probability  of  spare  availability  estimates  for 
Models  1,  2,  and  3.  Thus,  the  significance  of  the  difference  between 
the  probability  of  spare  availability  estimated  by  the  different 
simulation  models  can  be  determined.   In  addition,  the  significance  of 
any  difference  between  Model  1  and  renewal  theory  estimates  for 
probability  of  spare  availability  can  be  determined. 

Consider  a  random  variable  X  whose  density  function  is  given  by 

f(x|9).   The  likelihood  function  of  n  random  variables  I„  J.  ,„  I  , 

1  1   2      n 

denoted  by  L(X  ,  ...,  X  I 9) ,  has  the  joint  density  function  f   (...,  x  ) 
In1  J  '  x.       n 

(20).   The  likelihood  is  the  value  of  the  density  function,  which  for 

discrete  random  variables  is  a  probability. 

For  the  case  where  9  is  unknown  and  the  joint  density 

of  n  random  variables  is  f         (X, ,  ....  X  I 9) ,  the  maximum- 

xl    X  n1 

likelihood  estimator  of  9,  denoted  by  9,  is  the  value  of  9  which 

maximizes  the  likelihood  function  L(X  '  ,  ....  X  ' | 9) .   Thus,  for  many 

cases  the  maximum  likelihood  estimator  is  the  solution  to  (20) 

^-'   0  .  (2.118) 

Both  1.(9)  and  £nL(9)  have  a  maxima  at  the  same  value  of  9. 

The  probability  mass  function  for  the  number  of  failures  in  a 
random  sample  of  size  n  taken  from  a  binomial  distribution  is  given  by 

f(x|p)  =  px  q1_x  F{0>1}(x)  ,  (2.119) 

where 
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f{o,i}(x) 


0  <  x  _<    1 

otherwise 


0  <  p   <    1,    and 
q   =   1   -  p. 

The   likelihood   function  can  be  written  as 

n       x        1-x.  Zx       n-Zx 

L(p)    =nP1q        1=p1q  .  (2.120) 

1=1 

The  natural  log  of  the  above  equation  is 

Jin  [L(p)]  iji,  fln(p)  +  (n  -  Z  x.)  to(q)  .       (2.121) 

The  above  equation  differentiated  with  respect  to  p  gives 

3|jta[L(p)]  I        Zx.        n  -  Zx. 

.  (2.122) 


3p  p  q 

The  maximum  likelihood  estimator  p  is  obtained  by  setting  the  above 

expression  equal  to  0  and  solving  for  p.   Thus, 

Zx 
P  ■  —  =  x  .  (2.123) 

The  mean  of  x,  denoted  by  p  ,  is  found  by  taking  the  first  moment 

of  the  binomial  distribution  (21) 

n 
u  -  E  x  f(x)  ,  (2.124) 

x=0 

where  X  is  a  discrete  random  variable  and  f(x  )  is  its  density  function. 
The  mean  of  X  can  now  be  written  as 
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V    n   x/i   ii"X    r      n!       x,,   .n- 
"x  "  JQ  *[xj  P  C1"P}    "  J,  (x-DKn-x)!  P  (1"P> 


n 

I 
x=l 


-  np  I      \l:\\    p^U-p)"  .  (2.125) 


If  a  =  x  -  1  and  b  •  n  -  1 ,  then 


a=0 
where 


u  =np  I      M  pa(l-p)b-a  ,  (2.126) 


I      W  Pad-p)b-a  =  1.  (2.127) 

a=0  W 


Thus,  the  mean  of  X  is 

VK   =  np   .  (2.128) 

2 
The  variance  a  is  found  by  taking  the  second  moment  about  the 


origin  (21) 


,2=  I     x(x-l)  W  pX(l-p)n-x  (2.129) 

x=0        W 


Following  a  similar  procedure  as  outlined  for  the  mean  of  x,  we  find  the 
variance  is 

2     „   > 
a  =  np(l-p)  .  (2.130) 

The  maximum  likelihood  estimator  and  variance  can  now  be  used  to 
calculate  confidence  interval  estimates  for  the  probability  of  spare 
availability.   Let  p  be  the  proportion  of  simulations  when  the  spare  is 
available  at  some  time  t  .   If  repeated  samples  of  size  n  are  taken,  the 
number  of  simulations  for  which  the  spare  is  available  will  be  a  random 


59 


variable,  say  X.  For  fairly  large  values  of  n  the  distribution  of  X 

will  be  approximately  normal  (20).   If  the  distribution  of  X  is  normal, 

2 
we  say  the  random  variable  X  is  N(u,a  )  and 


X~  Ux 
Z ?T^^  N(0,1)  .  (2.131) 


As  was  previously  shown  u  =  np  and  a   =  /np(l-p) ,  therefore,  Eq.  (2.131) 
can  be  written  as 

X 

z  ■  *  -  'p  =  JLLL    .  (2.132) 

^np(l-p)    /p(l-p) 


The  standard  deviation  for  the  probability  of  spare  availability  at  t  , 

denoted  by  o(t  ) ,  can  be  estimated  by  setting  x  =  P    . ,  (x,  )  giving 
K  j  a  avail   k  °    s 


r*      ,i   (O  U-S       .,  (t,  )] 

avail   k'     avail   k  ' 


»(Tk)  "  /  "*"        S  g "^   -        ,  (2.133) 

where  M  is  the  number  of  simulations  run.   The  (1-ct)  confidence  interval 
is  given  by 


Prob  "  Za/2  <  S—%<   Za/2  "  l  '   °  •  <2-134' 


Pre,  -  Za/2  <  Jg-  <  Za/2  -!-..  (2.135) 

*■  M  J 


Thus,  the  100(l-a)%  confidence  interval  for  p  is  given  by 
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n  "  a/2  /   M    '  (2.136) 


Recall  that  the  maximum-likelihood  estimator  for  p  is  P     (t  ) 

avail  k 

therefore  the  100(l-a)%  confidence  interval  for  p  is  given  by 


P  CO   *  I  /Pavail(y[1-Pavail^ 

avail  lV  a/2      /  M  '  (2.137) 


where  Pava±1  (i"k)  is  the  average  probability  of  spare  availability 
determined  by  the  computer  code. 

2.6  SIMULATION  Computer  Code 

The  computer  code,  SIMULATION,  was  written  to  analyze  Models  1,  2, 
and  3  (see  Appendix  A).   The  computer  code  is  used  first  to  simulate  a 
random  number  over  the  interval  [0,1],  using  a  modified  version  of  the 
IMSL  subroutine  GGUBFS  (22) .   The  most  recent  repair/replacement 
interval  is  then  analyzed  to  determine  the  S.'s   by  computing  the  time 
difference  between  the  most  recent  failure  and  any  failures  prior  to 
that  using  Eq.  (2.79) . 

After  the  6.  values  have  been  computed,  the  value  of  &.'s 
distribution  function  is  computed  and  compared  to  the  random  number.   If 
the  value  of  the  random  number  is  greater  than  the  numerical  value  of 
the  distribution  function,  the  distribution  function  associated  with  6 
is  computed.   The  first  distribution  function  value  which  is  less  than 
or  equal  to  the  random  number  is  used  for  the  simulation  equation. 

The  simulation  equation  is  the  distribution  function,  where  the 
random  number  is  substituted  for  F(At  )  and  solved  for  At  ,  e.g.,  in 
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Eq.  (2.70).   Thus,  the  value  of  the  present  failure  time  interval  is 
computed. 

The  next  random  number  is  created  and  the  above  procedure  is 
repeated  until  the  sum  of  the  failure  time  intervals  exceeds  the  pool 
lifetime.   At  this  point,  a  "history"  of  simulated  plant  failures  has 
been  compiled.   The  failure  time  just  prior  to  the  failure  time  that 
exceeded  the  pool  lifetime  is  then  compared  to  the  pool  lifetime.   The 
difference  between  the  two  times  is  used  to  estimate  the  probability  of 
spare  availability  for  Model  1  and  fractional  downtime  for  Models  1  and 
3. 

In  addition,  the  last  time  the  spare  was  used  is  compared  to  the 
pool  lifetime  and  the  difference  is  used  to  estimate  the  probability  of 
spare  availability  and  fractional  downtime  for  Model  2.   Estimates  of 
the  probability  of  spare  availability  and  fractional  downtime  are 
obtained  by  averaging  the  probability  of  spare  availability  and 
fractional  downtime  over  a  large  number  of  histories  (e.g.,  for  this 
study  10,000  simulation  histories  were  considered  sufficient). 

2.7  Sample  Problem  For  a  Five  Component  System 

A  sample  problem  was  developed  to  compare  the  simulation  model 
results  to  renewal  theory  results.   The  first  version  of  the  example 
problem  considered  a  system  with  five  plants  in  the  spare  pool  and  a 
6  =  .02  (see  Fig.  2.13).   A  study  of  the  Fig.  (2.13)  indicates  that  the 
Model  1  results  for  probability  of  spare  availability  were  within  ±a  of 
the  renewal  theory  prediction  for  all  but  two  data  points.   In  addition, 
the  Model  3  results  were  not  different  from  the  Model  1  results.   Thus, 
for  small  values  of  S  the  numerical  solution  to  the  renewal  equation  can 
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Plant  component  failure  time. 

FIG.  2.13.  Probability  of  spare  availability  for  a  five  component 
system.  The  Model  1,  2,  and  3  simulation  results  are  compared  to 
the  renewal  theory  results  with  B  -  .02  and  10,000  simulations 
computed.  The  error  bars  indicate  ±a. 
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be  used  to  estimate  the  probability  of  spare  availability  for  the  Model 
1  and  Model  3  assumptions. 

The  second  variation  of  the  example  problem  considered  the  same 
system  as  before  but  by  using  simulation  Model  1.   The  number  of 
simulations  was  increased  by  a  factor  of  ten  (see  Fig.  2.14).   As 
expected,  the  Model  1  results  were  even  closer  to  renewal  theory 
prediction  than  the  previous  Model  1  results  using  fewer  simulations. 

The  third  variation  of  the  example  considered  the  same  system  with 
the  initial  seed  value  used  by  the  random  number  subroutine  changed  to 
check  for  the  "randomness"  of  the  simulated  failures  (see  Fig.  2.15). 
A  study  of  the  figure  indicates  that  the  fluctuations  in  the  probability 
of  spare  availability  occurred  at  different  failure  times.  This 
indicated  the  fluctuations  were  due  to  the  random  numbers  and  not  the 
logic  of  the  computer  code. 

In  the  fourth  and  final  variation  of  the  above  example  problem 

considered  a  system  with  five  plants  in  the  spare  pool  and  a  B  =  .A  (see 

Fig.  2.16)  was  considered.   Once  again,  the  Model  1  results  agreed  with 

the  renewal  theory  prediction.   The  Model  2  and  Model  3  results  were 

significantly  different  for  this  case.   With  the  decision  rule  used  in 

Model  2  (i.e.,  spare  is  used  only  if  immediately  available  at  the  time 

of  failure),  failures  can  occur  in  the  interval  (0,  T   ,)  and  the  spare 

repl  v 

would  still  be  available.   The  probability  of  spare  availability  was 
therefore  substantially  larger  for  times  greater  than  6  in  the  case  of 
Model  2  when  compared  to  Models  1  and  3.   Because  Model  3  assumes  a 
finite  spare  failure  rate,  the  probability  of  spare  availability  was 
different  for  Models  1  and  3. 
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FIG.  2.14.   Probability  of  spare  availability  for  a  five  component 
system.   The  Model  1  simulation  results  are  compared  to  the 
renewal  theory  results  with  B  =  .02  and  two  different  numbers  of 
simulations  run.   The  error  bars  indicate  ±o. 
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FIG.  2.13.  Probability  of  spare  availability  for  a  five  component 
system.  The  Model  1  simulation  results  are  compared  to  the 
renewal  theory  results  for  two  initial  seed  values  for  the 
random  number  generating  subroutine  with  10,000  simulations 
computed.  The  error  bars  indicate  ±o. 
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Model  1 
Model  2 
Model  3 
Renewal  theory 


2  3 

Plant  compongnt  failurg  timg. 

FIG:  2-1,i!;  P^b^ili-y  of  spare  availability  for  a  five  component 
system.  The  Model  1,  2,  and  3  simulation  re3ults  are  compared  to 
the  renewal  theory  results  with  8  =  .4  and  10,000  simulations 
computed • 
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The  simulation  models  and  renewal  theory  estimate  the  same 
probability  of  spare  availability  for  the  interval  (0,g).  This  result 
supports  the  hypothesis  that  the  spare-component  pool  (without  any 
simplifying  assumptions)  is  a  "pseudo-renewal  process."  That  is,  the 
probability  of  all  plants  in  the  system  operating  at  time  T  is  dependent 
upon  their  operational  status  in  the  interval  (t-B,t)  for  t  >  B. 
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Chapter  3 
ECONOMIC  MODEL  FOR  SPARE-COMPONENT  POOL 

3.1  Model  Formulation 

An  economic  model  that  accounts  exactly  for  all  possible  costs  and 
revenue  requirements  along  with  the  various  timing  of  these  cash  flows 
is  a  very  complex  problem.   A  simpler  model  utilizing  the  concept  of  a 
fixed  charge  rate  to  treat  capital  expenditures  was  used  to  analyze 
results  for  this  study.  To  evaluate  the  attractiveness  of  either 
management  plan  a  present  worth  cost  of  all  future  expenses  for  Plan  A 
and  Plan  B  are  estimated  and  then  compared.   The  plan  with  the  lower 
present  worth  cost  is  chosen  as  the  superior  option. 

3.1.1  Revenue  Requirements 

The  economic  model  was  based  on  the  concept  of  a  fixed  charge  rate 
FCR(n.m)  and  a  discount  rate  i,(n,m),  both  of  which  are  generally 
functions  of  time  (year  m)  and  the  plant  of  interest  (plant  n) .   In  this 
study  the  assumption  was  made  that  the  discount  rate  was  the  same  for 
all  plants  in  the  pool.   Thus,  the  discount  for  year  M  is  denoted  by 
Ijfrn).   Let  AR(n)  be  the  levelized  annual  revenue  requirement  over  m 
years  and  C(n,0)  be  an  expenditure  at  the  beginning  of  m  years.   The 
fixed  charge  rate  is  defined  such  that  FCR(n,m)C(n,0)  is  the  dollar 
amount  that  must  be  earned  by  the  end  of  each  year  over  a  period  of  m 
years  to  capitalize  the  expense  C(n,0)  for  plant  n  (23).   Thus,  AR(n)  is 
given  by 

AR(n)  =  C(n,0)FCR(n,m).  (3.1) 

The  FCR(n.m)  includes  considerations  for  return  on  investment,  return  of 
investment,  taxes,  depreciation,  etc.   Thus,  FCR(n.m)  gives  an  annual 
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levelized  capital  requirement  over  m  years  that  is  equivalent  to  the 
actual  requirements  obtained  by  a  detailed  analysis  of  the  many  factors 
influencing  utility  economics. 

The  total  revenue  required  R(n,0)  at  time  zero  to  cover  expense 
C(n,0)  without  any  further  revenue  is  the  present  worth  of  AR(n) ,  which 
is 


R(n,0)  =  AR(n) 


present 

worth 

factor 


(3.2) 


where 


present 

worth 

factor 


[H-id(D] 


[i+id(i)][i+id(2); 


+  [  1+i.  U)  ][  1+i.  (2)  ]  [  1  +  1.(3)]  + 


[i+id(D  ]  [i+id(2) ] . . .  [i+id(k-D  ]  [l+id(k)  ]Yk- 


Y,  is  the  fraction  of  the  last  year,  k,  over  which  the  revenue  must  be 


expenditure  that  occurs  somewhere  between  the  beginning  and  the  end  of 
year  k,  y^   will  be  a  fractional  value. 

The  notation  in  Eq.  (3.2)  can  be  simplified  to  give 


I  1  1 

R(n,0)  =  AR(n)  \    J   n  [1+i  (h)]"1 

J-l  h=l    d 


(3.3) 


The  results  of  Eq.  (3.1)  and  Eq.  (3.3)  are  combined  to  get  the 
relationship  between  initial  expenditure  and  initial  revenue  requirement 
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I   1  1 

R(n,0)  =  C(n,0)  FCR(n.k)   T   H  [l+i.(h)]"x  .    (3.4) 

i-1  h=l    d 

Therefore,  the  FCR(n.k)  divides  the  initial  expenditure  into  equal 
payments  over  k  years  such  that  the  sum  of  the  discounted  annual  revenue 
requirement  is  equal  to  the  initial  revenue  requirement. 

For  the  case  of  a  constant  discount  rate,  the  present  worth  of  all 
the  AR(n)  is 


I, 

R(n,0)  =  AR(n)  (P/A)  "  ,  (3.5) 

id 
where  (P/A).   is  the  present  worth  of  a  uniform  series  expressed  as 

i,    (l+!Jk-l 
(PM)k  "  " k  •  »•« 

Thus,  the  initial  revenue  requirement  for  an  initial  expenditure  of 
C(n,0)  is 


R(n,0)    =  C(n,0)    FCR(n,k)(P/A)    a   =  C(n,0)    FCR(n.k) 


<l+i/-l 


vi+vk 


(3.7) 


3.1.2  Revenue  Requirements  and  Total  Cost  Calculations  For  Management 
Plans 

In  this  section  a  method  is  derived  to  evaluate  the  cost  and 

revenue  requirements  for  the  two  component  management  plans .   The  two 

plans  are  compared  on  the  basis  of  the  present  worth  of  all  future 

expenditures.   Certain  types  of  expenditures  such  as  the  spare  component 

cost,  replacement  component  and  salvage  values,  are  evaluated  in  terms 

of  an  equivalent  revenue  requirement.   Other  expenditures  such  as 

operation  and  maintenance  costs,  replacement  energy,  and  reserve 

capacity  are  expensed  in  the  year  these  costs  are  incurred. 
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The  present  worth  of  costs  occurring  during  year  m  is  given  by  Eq. 
(3.4)  or  (3.7),  if  the  item  is  expensed.   In  contrast,  items  that  are 
capitalized  over  k  years  have  an  associated  FCR(k)  and  a  present  worth 
of  the  revenue  requirement.  For  an  item  that  is  capitalized  over  m 
years,  the  revenue  requirement  in  year  m  for  an  expenditure  C(m)  is 

R(n,m)  =  REV(n,k)C(n,m)  ,  (3.8) 

where  REV(n,k)  is  the  revenue-to-expense  conversion  factor  for  plant  n. 

The  present  worth  of  the  revenue  requlrment  at  time  zero  is  the 
revenue  requirement  in  year  k  discounted  over  m  years.   Thus, 

m  1 

R(n,0)  =  REV(n.k)  C(n,m)   II   [l+i,(j)]    =  REV(n.k)  C(n,0)  .   (3.9) 

j  =  l 

An  expression  for  the  revenue-to-expense  conversion  factor  is  found  by 
equating  Eq.  (3.4)  and  Eq.  (3.9)  to  give 

I   j  -1 

REV(n.k)  =  FCR(n.k)   [   n   [1+i  00]  .  (3.10) 

3-1  h=l 

The  revenue-to-expense  conversion  factor  for  plant  n  is,  therefore, 
independent  of  the  time  at  which  the  expense  is  incurred.   Thus,  the 
revenue  requirement  at  time  t'  is  given  by 

R(n.t')  =  REV(n.k)  C(n.t')  =  FCR(n.k)  C(n.t') 


r   j  -1 

I      I  [l+i.  (h)] 


^dv"/j   .  (3.11) 

J-l  h=l     a 


For  the  remainder  of  this  study  the  revenue-to-expense  conversion  factor 
dependence  on  the  capitalization  period  will  be  assumed  and  denoted 
simply  as  REV(n) . 
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3.1.3  Allocation  of  Plan  B  Spare  Component  Costs 

Many  of  the  costs  considered  in  this  project  are  shared  by  both 
Plan  A  and  Plan  B;  however,  costs  associated  with  the  spare  component 
are  specific  to  Plan  B.   These  spare  component  costs  include  initial 

purchase  cost  C    (0),  maintenance  and  storage  costs  C   ,  and  the  final 

apr  ■         ms 

salvage  value  C   .  Each  plant  must  pay  a  share  of  these  costs  based  on 
some  allocation  method. 

The  fraction  of  costs  assigned  to  the  n-th  plant  participating  in 
the  spare-component  pool  are  defined  by 

F(n,0)  =  fraction  of  initial  spare  component  purchase  cost  at  the 
beginning  of  the  pool, 

F(n,m)  =  fraction  of  storage  cost  in  year  m,  and 

F(n,T   )  =  fraction  of  spare-component  salvage  value  after  pool  has 

been  dissolved  after  T    years, 
max  J 

The  two  allocation  methods  considered  in  this  project  are  (1)  allocation 
by  protected  capacity  and  (2)  allocation  by  number  of  plants  in  the 
pool. 

Protected  capacity  method 

The  protected  capacity  allocation  method  assigns  a  value  to  the 
fraction  F  equal  to  the  ratio  of  the  capacity  protected  by  the  spare 
component  in  plant  n  to  the  total  capacity  protected  of  all  the  plants 
in  the  pool.   Thus  if  PMW(n)  denotes  the  protected  capacity  in  plant  n 
and  there  are  a  total  of  N  plants  in  the  pool,  the  allocation  fraction 
for  the  initial  component  purchase  cost  is 


«W0~Ftt-.  (3.12) 

I   PMW(j) 
j-1 
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The  allocation  fraction  for  storage  costs  would  only  be  assigned  to 
plants  that  were  operating  during  the  year;  or 

„,   s    PMW(n) 

F(n,m)  =  - —  ,  (3.13) 

I  PMW(j) 
j  =  l 

where  the  summation  is  over  only  the  K  plants  operating  during  year  m. 

For  this  method,  the  salvage  fraction  F(n,T   )  is  taken  equal  to 

max 

F(n,0). 

Number  of  plants  method 

A  simpler  allocation  method  and  the  one  used  in  the  later  examples 
is  to  assign  costs  on  the  basis  of  the  number  of  plants  in  the  pool. 
That  is 

F(n.O).  =  F(n,Tmax)  =  I  ,  (3.14) 

where  N  is  the  number  of  plants  in  the  pool  and 

F(n,m)  =  irW   *    0<N'(m)<N,         (3.15) 

where  N' (m)  is  the  number  of  plants  in  the  pool  that  operate  during  year 
m.   Regardless  of  which  allocation  method  is  used,  no  annual  charge  for 
pool  expenses  is  assigned  to  non-operating  plants.   In  this  manner 
plants  could  be  allowed  to  enter  the  pool  after  its  formation  as  well  as 
leave  the  pool  before  the  pool's  dissolution. 

3.2  Classification  of  Cost  Components 

The  types  of  costs  associated  with  both  management  plans  can  be 
divided  into  four  separate  classifications:   1)  failure-dependent 
variable  costs,  2)  plant-operation  variable  costs,  3)  failure-dependent 
fixed  costs,  and  4)  annual  fixed  costs  (see  Table  3.1). 
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Table  3.1  The  four  classifications  of  costs  associated  with  both 
management  plans  and  the  costs  in  each  category. 


Types  of  Costs 


Plan  A 


Plan  B 


Failure-dependent  variable  costs 
-replacement  energy 

Plant-operation  variable  costs 
-operation  and  maintenance 

Failure-dependent  fixed  costs 

-component  repair  or  replacement 
-used  component  salvage  value 
-failed  component  salvage  value 
-spare  used  as  temporary  substitute 

Annual  fixed  costs 

-reserve  capacity 

-spare  purchase 

-spare  storage  and  maintenance 

-spare  salvage  value 


X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

The  following  sections  describe  methods  of  estimating  these  costs 
for  each  classification.   Costs  classified  as  types  1,  2  and  3  are 
differential  costs  incurred  in  the  repair/replacement  interval  and  they 
are  present  valued  to  the  time  of  failure,  t'.   These  differential  costs 
are  then  summed  over  the  interval  (0,t)  to  estimate  the  present  worth  of 
the  type  costs.   In  contrast,  type  4  costs  are  incurred  at  the  end  of 
each  year  and  represent  the  value  at  that  time.   These  costs  are  then 
summed  over  the  interval  (0,k)  to  estimate  the  present  worth  of  that 
type  cost,  where  k  is  an  integer  value  of  years  such  that  k  -  1  <  t  <   k 
(e.g.,  if  t  =  2.3,  then  k  =  3) .   Thus,  the  present  worth  costs  are  an 
estimate  of  the  total  costs  for  plant  n  over  the  interval  (0,t).   In 
general,  the  value  of  t  will  be  equal  to  the  spare-pool  lifetime. 
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3.3  Total  Plan  Costs 

The  total  costs  under  the  two  managerial  plans  can  now  be 
estimated.  The  present  worth  of  all  expected  costs  at  time  t  and  at 
plant  n  under  plan  I   is 

PWTOTAL^n.t)  =  PWFVC^n, t)+PWOM£ (a, O+PWFFC^n.O+PWAFC  (n,t)  ,  (3.16) 

where  PWFVC  is  the  present  worth  of  the  failure-dependent  variable 
costs,  PWOM  is  the  present  worth  of  the  plant-operation  variable  costs, 
PWFFC  is  the  present  worth  of  the  failure-dependent  fixed  costs,  and 
PWAFC  is  the  present  worth  of  the  annual  fixed  costs. 

3.4  Failure-dependent  Variable  Costs 

The  variable  costs  at  the  time  of  failure  are  dependent  on  the  time 
of  failure  and  the  length  of  each  shutdown  period.   In  this  section  the 
general  case  is  analyzed,  followed  by  several  special  case  solutions  for 
the  different  management  plan  total  variable  costs.   In  general, 
numerical  techniques  must  be  used  to  evaluate  the  integrals  in  the 
expressions  for  these  costs. 

The  total  variable  costs  of  a  failure  are  dependent  on  the 
management  plan  the  plant  is  operating  under.   The  Plan  A  costs  are  the 
present  worth  of  the  variable  costs  that  are  paid  by  the  k-th  plant  over 
the  interval  (0,t),  and  in  general,  are  a  function  of  the  downtime 
following  a  component  failure.   The  Plan  B  costs  are  the  present  worth 
of  the  variable  costs  that  would  be  paid  by  the  k-th  plant  over  the 
interval  (0,t)  if  the  plant  were  a  member  of  the  spare  pool.   These 
costs  are  the  weighted  average  of  the  variable  costs  for  times  when  the 
spare  is  not  available  (Plan  A  costs)  and  for  times  the  spare  is 
available  and  used.   In  general,  the  Plan  B  costs  incorporate  some 


76 

fractional  portion  of  the  Plan  A  costs  and  are,  therefore,  also  a 
function  of  the  downtime  as  well  as  the  probability  of  spare 
availability  at  the  time  of  failure. 

3.4.1  General  case 

In  this  section  a  model  is  derived  for  estimating  the  Plan  A  and 
Plan  B  failure-dependent  variable  costs.   The  present  worth  of 
failure-dependent  variable  costs  at  plant  n  due  to  a  failure  at  t1  can 
be  expressed  as 


n-th  plant  fails 
fvc(n,t')  =  Prom  in  unit  time 
about  t' 


n-th  plant  is 
operating 
at  fc' 


fraction  of  the 
'"\  repair /replacement 
time  the  plant  is 
shutdown  or  derated 


cost  of 

■ 
present 

failure 

worth 

at  t ' 

discount 
factor 
L        ) 

Recall  from  renewal  theory  that  g(t')  Is  defined  as 


any  plant  is  ' 
ProM  operating  (without 
the  spare)  at  t' 


thus 


=  g(t'), 


(3.17) 


a  plant  failure 
Prob-j  in  unit  time 
about  t ' 


n-th  plant 

is  operating \   =  Prob 

at  t' 


a  plant 
failure  in 
unit  time 
about  t' 


n-th  plant 
•  ProM  is  operating  > 
at  t' 


*  g(t'). 


(3.18) 


Define  C(n.t')  as  the  dollar  rate  per  unit  downtime  at  plant  n  and  time 
t '  making 


cost  of  a 
failure  at  t ' 
in  plant  n 


=  C(n,t')  T 


repl' 


(3.19) 
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In  general,  the  cost  per  unit  downtime  will  increase  over  time  by  so 
yearly  escalation  rate  making 


C(n.t')  =  C( 


C-l 
n 


[1  +  I(j)]   [1  +  I(m)]  m, 


(3.20) 


where  m  =  INT(t')  (i.e.,  m  is  the  largest  integer  value  <  t'),  y  ■ 

~       m 

t'  -  m,  and  I(m)  is  the  escalation  rate  in  year  m.  The  present  worth 
discount  factor  converts  a  future  dollar  amount  into  an  equivalent 
present  dollar  amount.  Let  F(t')  be  the  future  value  of  P(0)  present 
dollars  making 


P(0)  =  F(t') 


m-1 

n  [i+i.(j)]" 


[l+id(m)] 


(3.21) 


where  m  =  INT(t'),  y^   =  t'-m,  and  id  is  the  discount  rate  in  year 
For  the  case  of  a  constant. discount  rate  P(0)  becomes 


P(0)  =  F(t')  (1+i.) 
d 


■t ' 


(3.22) 


The  present  worth  of  the  failure-dependent  variable  costs  due  to 
failure  at  t'  can  now  be  written  as 


fvc(n.t')  =  \   g(t')  fDT(t')  C(n,0)  Tre  x   I(t') 


(3.23) 


where  fDT(t')  is  defined  as  the  "fractional  downtime"  defined  as  the 
effective  fraction  of  the  repair/replacement  time  that  the  plant  is 
shutdown  or  derated  while  recovering  from  a  component  failure.   The  term 
I(t')  is  given  by 


I(t')  = 


V  imm 


j-i  ^i+id^ 


1+I(m) 
l+id(m) 


(3.24) 
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where  m  =  INT(t')  and  y  =  t'-m,  and  represents  the  present  worth 
factor. 

The  total  present  worth  of  the  failure  dependent  variable  costs  is 
the  sum  of  all  the  differential  costs  at  t',  or 

t  t 

FWFVC(n.t)  =  J  dt'  fvc(n,t')  =  j  dt'X  g(t')  ^(t')  C(n,0)  T^  ±   I(t') 
0  0 

t 
=  \   C(n,0)  Trepl  J  dt'  g(t')  fDT(t')  I(t').  (3.25) 

For  the  case  of  constant  escalation  and  discount  rates,  this  result 

simplifies  to 

t  t' 

PWFVC(n.t)  =  A  C(n.O)  T^  J  df  g(t')  E^Ct')  [|§J   .    (3.26) 

0  d 

For  Plan  A,  *j,T(t')  =  1,  so  that 

t 
PWFVCA(n,t)  =  X   C(n,0)  Trepl  J  dt'  g(t')  I(t')  ,       (3.27) 

0 

while  for  Plan  B,  fDT(t')  is  some  function  of  time  t'.   In  general,  Eqs. 
(3.25)  or  (3.26)  must  be  evaluated  numerically  because  of  the  piecewise 
construction  of  g(t'),  the  variable  discount  and  escalation  rates,  and 
the  fractional  downtime.  Only  for  very  simple  situations  in  which 
analytical  representations  of  the  quantities  are  available  is  it 
possible  to  evaluate  the  integrals  in  Eqs.  (3.27)  or  (3.28) 
analytically.   In  the  following  sections,  PWFVC  (n,t)  is  evaluated  for 
two  different  system  operating  conditions  and  possible  forms  of  g(t'). 


79 


3.4.2  Case  With  Constant  Plant  Component  Availability  And  Fraction 
Downtime 

In  this  section  the  failure-dependent  variable  costs  for  Plan  A  and 
Plan  B  are  evaluated  under  the  following  assumptions: 

1)  The  annual  escalation  and  discount  rates  (I  and  i.)  are 

d 

constant. 

2)  The  probability  of  spare  availability  and  the  probability  of  a 
plant  operating  without  the  spare  at  a  given  instant  in  time 
are  both  determined  by  the  asymptotic  renewal  theory  result. 

3)  The  spare-component  failure  rate  is  zero. 

4)  The  fractional  downtime  is  a  constant. 
Under  the  above  conditions,  g(t')  becomes 


g(t')  ■  x«l (3-28> 

repl 


so  that  the  present  worth  of  the  Plan  A  variable  costs  are 

fvcA(n,f)    =   X  df    ['^ ]    a£'    C(n.O)    T  ,  (3.29) 

v  repl-1  r 


where   a  =    (1  +  1) /(1+iJ  . 
d 

The  present  worth  of  the  Plan  A  variable  costs  for  plant  n  over  the 
interval  (0,t)  becomes 

t  t 


,  at  C(n,0)  T   ,.  (3.30) 
Trepl-1  rePX 


PWFVCA(n(t)  -  ^dffvc^n.t')  -  fdfA^ 


Evaluation  of  the  integral  gives 


C(n,0)  T     X    (a'-l) 

PWFVC.(n.t)  =  ,  .  ^P1 ■ (3  in 

Av         to(a)(l  +  A  T   ,)     -  U.Ji; 

repl 


80 

The  present  worth  of  the  Plan  B  variable  costs  for  plant  n  over  the  same 
Interval  (0,t)  Is 


PWFVCB  (n,t)  =  J  dt'fvcB(t') 


=  dt'X 
0 


1+XT 


repl 


a   C(n,0)  T     fnT(t'),    (3.32) 


repl  DT 


where  a  =  (l+I)/(l+i  ).  This  result  reduces  to 


C(n,0)  I     X  (aC-l) 

'"S^  ■   Ma)   !'+  XT   ,)  fDT" 

repl 


(3.33) 


3.4.3  Case  with  Constant  Fraction  Downtime 

Now  consider  the  system  described  by  the  following  conditions: 

1)  The  annual  escalation  and  discount  rates  are  constant. 

2)  The  probability  of  a  plant  operating  without  the  spare  at  a 
given  instant  in  time  is  given  by  g(t')  where 

-Xt' 


g(t')  = 


1 


1+X  T 


repl 


0  <  t'  <  T 

—  repl 


T   ,  <  t  <  »  . 
repl     — 


(3.34) 


3)   The  fractional  downtime  is  a  constant  value  equal  to  the 
asymptotic  renewal  theory  result. 
The  present  worth  of  the  Plan  A  variable  costs  is 


fvc  (n,t')  = 


X  e'^V  C(n,0) 


repl 


,   0  <  t'  <  T 

—  repl 


1+XT 


repl 


t' 
a  C(n,0)  T   ,   ,   T   ,  <  t'  < 
repl      repl      - 


(3.35) 
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where  a  =  (l+I)/(l+id) .  The  present  worth  of  the  Plan  A  variable  costs 
for  plant  n  over  the  interval  (0,t)  is 


dt'  X  e"Xt'at'c(n)0)  T   ,     ,    0<t<T   , 
)  repl  -  repl 


PWFVCA(n,t)  = 


dt'  X  e"Xt'at'c(n,0)  T 


(3.36) 


repl 


+  Jo"'  X   M"^  at'C(n,0)  T"p^  '  W*" 


The  above  equation  is  solved  using  integration  by  parts  resulting  in 


PWFVCA(n,t)  =  ■ 


C(n,0)  Tre   X  (»t*"it-I)  e"Xt 


X  +  (in(a) 

C(n,0)  T   .  X  (aVU-l)  e"Xt 
repl 


,   0  <  t  <  T 


repl 


(3.37) 


X  +  Jin  (a) 


C(n)  Trgpl  X  (.*-!) 

2n(a)(l+X  T   ,) 
repl 


■   T     <  t  <  »  . 
repl     - 


For  a  constant  fractional  downtime  the  present  worth  of  the  Plan  B 
variable  costs  for  plant  n  over  the  interval  (0,t)  is 


making 


PWF?CB(n,t)  =  fDT  PWFVCA(n,t) 


C(n,0)  Trepl  X    UVXt-l) 


(3.38) 


PWFVC„(n,t) 

D 


X   +  Jin  (a) 


C(n'0)  Trepl  A  (ate"Xt-l) 
X  +  £n(a) 


C(n,0)  Trepl  A  (*  -» 


f    ,   0  <  t  <  T   , 
DT         -  repl 


(3.39) 


f„T  •    T 1   <  *  <  •  • 


teU)  (l+X  Tre  x)      LDT  '    repl 
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3.4.4  Replacement  Energy  Costs 

One  important  cost  incurred  when  a  plant  is  forced  to  be  derated  is 

the  replacement  energy  costs.   In  general,  the  spare  will  not  have  the 

same  rated  output  as  the  failed  component  it  is  replacing.   If  the  spare 

is  rated  higher  than  the  failed  component,  no  replacement  energy  costs 

are  incurred  during  the  time  in  which  the  spare  is  installed.   In 

contrast,  if  the  spare  is  rated  lower  than  the  failed  component  the 

plant's  output  capacity  is  decreased  from  MW    to  [MW    -  MW    1. 

max       max     down 

where  MWjown  is  the  output  capacity  lost  due  to  the  component  failure. 

In  the  case  of  Plan  A  the  decrease  in  output  capacity  is  MW 

max 

Replacement  power  must  be  purchased  from  outside  the  plant  to  make 

up  for  the  loss  in  output  capacity  as  a  result  of  the  component  failure. 

The  system  replacement  energy  costs  C   (n,t)  in  the  interval  (t,  ,  t,  + 

en  k   k 


repl)  ls  8iven  by 


Cen(n,t)  "  fc«p<n,t)  ""down*"*  W»'t}  (3-40) 


where  the  three  terms  are  defined  by: 


fc   (n,t)  =  the  plant  capacity  factor  averaged  over  the  interval 

(\>    tk  +  Trepi)  or  the  ratio  of  the  scheduled  plant 

output  (prior  to  the  failure)  to  the  maximum  possible 

plant  output  for  time  T 

repl 

MWdown^n^  =  the  Plant  output  which  is  lost  because  of  the  failed 

component.   If  the  failed  component  is  necessary  for 

the  plant  to  operate  at  any  output  then  MW,    is  the 

down 

plant  rated  capacity.   In  general,  MW,    is  the  failed 

down 

component  capacity. 
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?   (n,t)  =  the  differential  cost  of  replacement  energy  ($/down 
time)  in  the  interval  (t  ,  t.  +  I    ).   The  average 
cost  of  replacement  energy  is  the  difference  between 
the  cost  of  purchasing  power  contracted  from  some 
outside  source  and  the  cost  of  producing  power  at  the 
effected  plant.   An  allowance  is  made  for  escalation  by 
letting 

rm-1  -i  Y 

Cen(n,t)  =  Cen(n,0)  I  "  [1+Ien(J) ]J [1+Ien(m) ]   '   (3>41) 

where  I   (m)  is  the  escalation  rate  in  year  m,  m=INT(t) 
en 

and  Y  =  t  -  m. 
m 

For  Plan  A,  the  downtime  is  equal  to  the  new  component  repair/ 

replacement  time  and  MW     is  equal  to  the  full  plant  capacity,  MW(n) . 

For  Plan  B,  the  downtime  is  equal  to  MW.    which  is  defined  by 

down  ' 

Max(0,MW(n)  -  MW   )  where  MW    is  the  spare  capacity, 
spr         spr        r  r         ' 

If  the  results  of  Eqs.  (3.37)  and  (3.40)  are  combined,  the  present 
worth  of  the  failure-dependent  variable  costs  under  Plan  A  for  plant  n 
over  the  interval  (0,t)  for  the  general  case  are 


PWFVC  (n,t)  =  T     MW(n)  S   (n,0)  X       dt'  g(t')  f    (n,t')  I   (t'),  (3.42) 
a         repj.        en        J  ..  cap        en 


where 


en 


m-1  [1+1   (j)] 
en  J 


j=i  ri+VJ>! 


l+I   (m) 
en 


t  l+ld(«) 


and  m  =  INT(t)  such  that  y  ■  t'  -  m. 
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The  above  result  can  be  simplified  considerably  by  assuming  a 
constant  discount  rate,  constant  escalation  rate,  and  a  constant  capacity 
factor.  For  this  case 

PWFVCA(n,t)  =  Trepl  fcap(n)MW(nKen(n,0)X  j  df  g(f  )  UW   .  (3.43) 


If  the  spare  is  not  rated  for  full  capacity  output  of  the  plant, 
the  utility  must  purchase  energy  during  the  time  the  spare  is  installed. 
The  models  presented  in  Chapter  2  are  used  to  determine  the  fraction  of 
the  repair/replacement  interval  for  which  replacement  power  must  be 
purchased.  Thus,  the  present  worth  of  the  Plan  B  total  variable  costs 
at  the  time  of  failure  can  be  estimated  by 


PWVCB(n,t)  =  Trepl  Cen(n,0H   MW(n)A  j^df  g(t '  )f  cap(n,  t  '^(t ')  Iffi(t ' ) 


""down  ALdt'  •<«*>  «„,<«*'>  [1-fDT(t,)J  lm^">\ 


(3.44) 


where 


m-1 
J-l 


1+Ien^ 


l+idU) 


1+1   (m) 
en 


l+id(m) 


and  m  =  INT(t)  such  that  v  -  f  -  m. 
m 

Recall  MW     =  MW(n)  -  MW   ;  thus,  the  Eq.  (3.44)  can  be  written  as 
down  spr 

PWFVC  (n,t)  =  T     r   (n,0)X  I   MW     J  dt'  g(t')f    (n,t')I   (t 
repj.  en  down  J  n  cap       en 


spr 


df  g(t')  f   (n.f )  f-T(f )  I   (f ) 


cap 


DT' 


(3.45) 
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Equation  (3.45)  can  be  simplified  further  by  assuming  a  constant  capacity 
factor,  discount  rate,  and  escalation  rate  so  as  to  give 

PWFVCB(n,t)  =  Trepl  ?en(n,0)X  f^n)   MW^  \  it<    l(f  >[^* 

+  MWsprJ0dt'  8^>  fDT(£')[lTifJ     •  C3.46) 


Equations  (3.45)  and  (3.48)  can  be  solved  analytically  in  the  case  of 
constant  g(t')  and  f   (t').  The  results  are 

POTVCA(n,t)  =  <(n)  i|_=Ii  MW(n)  (3.47) 


and 


PWVCB(n,t)  -  <(n)  (^~^  [MW(n)  -  (l-Cggfe))  H»_J, 


where  <(n)  =  T     f   (n)  ?   (n,0)  X  g,  and  a 


spr 

1+1 
en 


repl  capv"'  "en"""'    "  s«  ■""  *    TTiT 

d 


3.5  Plant-operation  Variable  Costs 

Expenditures  during  the  time  of  plant  operation  include  operation, 
maintenance,  testing,  and  insurance  costs.  For  the  purposes  of  this 
study  it  is  assumed  these  costs  are  the  same  (on  a  cost  per  unit 
operation  time  basis)  for  the  original  component,  repair/replacement 
component,  and  spare  component  installed  at  a  particular  plant.   The 
variable  costs,  denoted  by  vcom,  are  a  function  of  the  particular  plant 
and  the  operation  time. 

The  present  worth  of  the  variable  cost  during  the  time  of  plant 
operation  in  unit  time  about  t '  is  expressed  by 


f        If  )    present  \  r 

plant  is    cost  per  unit   worth  plant  fails 

vcom(n,t')  =  ProM operating W operating    W discount [■  +  Prob  \ln   unit  time 
at  t1       time  factor  J        about  dt' 
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plant  is 
operating \ 

at  t' 


me 

is   H_ 


I  fraction  of  downtime 
.v    ^,  u  .   I  J  downtime 

during  which  spare  is  t*\ 
.-    .  .     per  failure 

installed  and  operating   ' 


Iroresent  "\ 
worth 
•i discount ? 
factor 


The  present  worth  of  the  Plan  A  plant-operation  variable  costs  for  plant 
n  in  unit  time  about  t*  can  now  be  written  as 


vcom.(n,t')  =  g(t')  ;   (n,0)  I   (t?) 
a  om       om 


where 


^ Ct')  = 

om 


-i  [l+i   a)] 

T      om 


j=i  [!+*,«>] 


l+I   (m) 
om 


l+id(m) 


(3.49) 


and  m  =  INT(t')  such  that  Y  -  t'  -  m.   The  present  worth  of  all  Plan  A 

m 

plant-operation  variable  costs  for  plant  n  over  the  interval  (0,t)  is 


PWOM  (n,t)  =   dt'  vcom  (n.t1)  =  C   (n,0)    dt'  g(t')  I   (t?)  .   (3.50) 


The  present  worth  of  the  Plan  B  plant-operation  variable  cost  for 
plant  n  in  unit  time  about  t'  is 

vconL(n,t')  =  g(t')  ?   (n,0)  1   (t') 
&  om       om 

♦  X  g(t')  [l-fDT(f)]  T    ^(n.O)  I  («.)  .    (3.51) 


The  present  worth  of  all  Plan  B  plant-operation  variable  costs  for  plant 
n  over  the  interval  (0,t)  is 

f  t 


PWOM^n.t)  =   dt'  vcomB(n,t')  =  C 


(n,0) 
en 


dt-  g(t')  i   (f ) 

om 


A  Trepl  jdt'  •<*')  U-^f)]  iom(f) 
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Com(n,0)   (WI    >   dt'  g(t')  I   (f) 


-XTreplJ0dt's(f)  fDT(f)Iom(f) 


(3.52) 


For  the  special  case  of  a  constant  g(t),  discount  rate,  escalation 
rate,  and  fractional  downtime,  the  costs  for  the  two  plans  become 


PWOMA(n,t)  =  Com(n,0)  gf  « ■  ^<f)  -  Com(n,0)  g(t)  (f^ 


(3.53) 


and 


PWC*L(n,t)  =  5   (n,0)  g  Ul+X  T   ,  - 

om  repl    DT   repl'  )  Q" 


W  *   .)   dt'  I  CO 


=  ?   (n,0)  g  (1  +  X  T   ,  -  f  X  T   ,)  f-^-^L 
°m  repl    DT   repr [in (a) 


(3.54) 


where  a  =  (1+1  )/(l+i,), 
om      d 


3.6  Failure-Dependent  Fixed  Costs 

Fixed  costs  occur  at  the  time  of  a  component  failure  and  are 

independent  of  the  length  of  time  the  plant  is  subsequently  shut  down  or 

derated.   These  costs  include  the  cost  of  repairing  the  failed  component 

or  purchasing  a  replacement  component,  C    (t) ,  and  the  salvage  value 

comp  ° 

of  the  failed  component  C„   (t) . 
f  sv 

Plan  B  costs  also  include  the  cost  of  taking  the  spare  component  to 
the  failed  plant,  installing  it,  and  returning  it  to  the  pool  when  the 


88 


repaired /replacement  component  arrives,  denoted  by  C    (t).  A  detailed 

temp 

description  of  the  failure-dependent  fixed  costs  is  provided  below. 

3.6.1  Repaired/replacement  Component  Costs 

The  repair/replacement  cost  at  time  t  for  plant  n  includes  the 

total  repair/replacement  and  installation  costs  for  the  component  such 

as  capital  cost,  repair  cost,  shipping,  removal  of  the  failed  component, 

etc.   The  utility  has  the  option  to  repair  the  failed  component  instead 

of  purchasing  a  new  component  in  which  case  C     is  the  total 

comp 

restoration  costs  associated  with  returning  the  failed  component  to  the 
"good  as  new"  condition.   An  allowance  is  made  for  the  escalation  in  the 
repair/replacement  cost  by 


m-1 

C (n.t')  =  C    (n,0)  \    n  [1+1    (j)] 

comp         comp       •  =  !    comp 


[1  +  t (m)]  m  ,  (3.55) 

comp 


I 

where  I    (m)  is  the  escalation  rate  in  year  m,  m  =  INT(t'),  and 

Y  ■  t'-m. 

m 

The  repair/replacement  cost  is  assumed  to  be  identical  for  both 

Plans  A  and  B.   In  the  event  the  spare's  rated  capacity  is  identical  to 

that  of  the  failed  component,  the  affected  plant  may  elect  to  keep  the 

spare  and  place  the  repair/replacement  component  in  the  spare  pool.   As 

a  result,  the  transportation  costs  may  be  different  for  Plan  A  and  Plan 

B  and  a  distinction  must  be  made  between  the  two  plans  for  C 

comp 

In  this  study  the  purchase  cost  was  assumed  to  be  capitalized  by 
the  utility.   Thus,  C     is  multiplied  by  the  revenue-to-expense 
conversion  factor  defined  in  Section  3.1.   The  value  of  REV(n)  is 
determined  for  the  affected  plant  and  the  year  in  which  the  purchase 
cost  occurs. 
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3.6.2  Failed  Component  Salvage  Value 

The  salvage  value  of  the  failed  component  is  zero  in  the  case  that 
the  failed  component  is  repaired.   An  allowance  is  made  for  the 
escalation  in  the  failed  component  salvage  value  by 

Cfsv(n't?)  "  Cfsv(n'0)  F^W*"]  t1+Ifsv(m))]Y,°  '    (3-56) 

where  I,   (m)  is  the  escalation  rate  in  year  m,  m  =  INT(t'),  m  and  y     = 
xsv  m 

t'  -  m.   Similar  to  the  repair/replacement  cost,  the  failed  component 
salvage  value  is  a  capitalized  expense. 

3.6.3  Used  Plant  Component  Salvage  Value 

A  used  but  operable  plant  component  also  has  some  salvage  value  at 
the  end  of  the  pool.   In  general,  the  salvage  value  will  depend  on  the 
size  and  age  of  the  component.   Like  the  component  repair/  replacement 
cost  and  the  failed  salvage  value,  the  used  salvage  value  is  capitalized 
by  using  a  revenue-to-expense  factor.   An  allowance  is  also  made  for  the 
escalation  in  the  salvage  value  of  the  used  component  by 


where  I    (»)  is  the  escalation  rate  in  year  m,  m  ■  INT(t'),  and  y 
pcsv  u 

t'  -  m. 


3.6.4  Cost  of  Installing  a  Spare  Temporarily 

For  Plan  B  only,  there  is  an  additional  cost  associated  with 
transporting  the  spare  component  to  the  failed  plant,  installing  it  and 
transporting  it  back  to  the  pool  storage  location.   If  the  spare 
capacity  is  identical  to  the  failed  component,  then  it  is  assumed  that 
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the  spare  remains  installed  in  the  plant  and  the  spare  transportation 
cost  is  neglected.   Like  the  other  fixed  costs,  the  temporary 
installation  costs  escalate  as 


W^^    =   Wn'°>     (^W^1)     [1+Wm)lYm'       »-58> 


where  I (m)  is  the  escalation  rate  in  year  m,  m  =  INT(t')  and  v 

temp  i  [J, 

t*  -  m. 


3.6.5  Present  Worth  of  Failure-Dependent  Fixed  Costs 

The  present  worth  at  time  t'  of  the  failure  costs  for  Plan  A  can  be 
expressed  as 

plant  n  fails 


ffc  (n,t')  =  Prob  ^in  unit  time 
A 

about   t' 


n-th  plant 
is  operating 
at   t' 


II    (present  \    f  } 

cost  of  a       worth  revenue-to- 

failure      rediscount  [""jexpense  con-      f 
at   t1  factor  version  factor 


where 


ffcA(n,t')  =  A  g(t')  REV(n) 


m-1  rl+l    (j) 

i   (t»)  =  n  comP 

compU  >         I  ^  (  l+id(j) 


1+1    (m) 
comp 

l+id(m) 


(3.59) 


If   (f)  = 
fsv 


"J1  f1+Ifsv(J 


J-1  1 1+id^ 


1  +  1,:   (m)-,Ym 
fsv 

l+id(m) 


m  =  INT(t'),  and,  f     =   t'-  m.   The  present  worth  over  the  interval  (0,t) 
of  the  failure  dependent  fixed  costs  include  the  summation  of  all  the 
failure  costs  less  the  revenue  from  the  sale  of  the  used  component. 
Thus,  the  present  worth  of  Plan  A  failure-dependent  fixed  costs  for 
plant  n  over  the  interval  (0,t)  is 
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PWFFC  (n,t)  -    df  ffc  ( 

1  Jo     A 


)f  "1  (  (/present  ■ 

component  is  1  revenue  from  worth 
-  Prob  <  operating   Msale  of  used U discount 
at  t         component     factor 


componen 


0 


PWFFC    (n,t)    -  dt'    ffc.(n.t')      -  g(t)    C       „(n,0)    REV(n)    I  (t) 


pcsv 


pcsv 


comp 


REV(n)    \X  (n,0)        dt'    g(t')    I  (t') 


comp 


*fsv' 


"   *   Sf      (".<»        d<'    8(f)    t.      (f) 


fsv 


where 


g(t)    t  (n,0)    I  (tU    , 

pcsv  pcsv 


I      (t) 

pcsv 


m-1    rl+I  (j) 

IL      I        Pcsv 


jll    I   1+id<J) 


1+1  (m)/ra 

pcsv 

1+ijGn) 


(3.60) 


The  present  worth  at  time  t '  of  the  failure  costs  for  Plan  B  are 
expressed  in  the  same  manner  as  for  Plan  A  with  an  additional  term  that 
accounts  for  the  temporary  installation  of  the  spare  as 


ffc  (n,t')  =  X  g  \z (n,0)  REV(n)  I    (t') 

1  c  omp  c  omp 


Cf8v(n.f)  REV(n)  I^Ct') 


P(n)  ?fomn(n,0)  I    (t*)| 
temp       temp    J 


where 


I    (f )  = 
temp 


m-1  fl+I    (j) 

n  temP 


l+I    (m)iY™ 
temp 

l+id(m) 


(3.61) 
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m  =  INI(t'),  and  y     =  t1  -  m.   The  term  p(n)  is  defined  by 


P(n) 


1 1   MW(n)  -  MWspR  f   0 
0   MW(n)  -  MWspR  =  0 


Thus,  the  present  worth  of  the  Plan  B  failure-dependent  fixed  costs  over 
the  interval  (0,t)  is 

It 
X?    (n,0)    dt'  g(t')  I    (t') 
comp      j  •■>  »   comp 


X?fsv(n,0)  j^dt'  g(t')  Ifsv(f) 


g(t)  C     (n,0)  I    (t) 
pcsv      pcsv 


X   P(n)  StemP(n'0)  Ldt"  8<f)  Itemp(f )  .  (3.62) 


For  the  special  case  of  constant  operating  probability  (i.e.,  a 
constant  g(t)),  escalation  rates,  and  discount  rates  the  analytical 
evaluation  of  Eqs.  (3.62)  and  (3.64)  yields 

PWFFV-f    ■  «!^j    "  <2(&)    "  <3    (&)  »•«) 


and 


t 


a  -1    )        .      bu-l  „    fcc-l   I     ,  fdC-l    ] 


"■**<*•«    "  <llin(aTj    "  <2[fnWJ    "  <3[Ut}    +  <4    ^J    •        (3.64) 


where 


Kl  =  X   ScomP(n'0)   REV(n)   «• 


K2   =  X   Cfsv(n,0)    REV(n)    8. 
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K3  =  Scsv(n'0)  REV(n)  «■ 


K4  =  p(n)  ?temp(n'0)  8' 


a  = 


b 


l+r 

comp 

1+id 


1  +  IF 

ISV 


11,  |   ' 


1  +  1 

pcsv 

1+1d 


temp 
1+1d 


,  and 


3.7  Annual  Fixed  Costs 

Certain  fixed  costs  are  incurred  by  the  generating  station  even  if 
no  failures  occur.   Such  fixed  costs  are  considered  year-end  expenses 
and  must  be  estimated  in  annual  increments.   The  assumption  is  made 
that  the  plants  will  contract  for  reserve  margin  capacity  and  pool 
membership  on  a  yearly  basis.   Unlike  the  energy  replacement  cost,  the 
plant  reserve  capacity  cost  is  contracted  and  paid  for  on  a  yearly 
contract  basis.  Thus,  the  system  reserve  cost  is  a  discrete  function  in 
time  and  can  be  evaluated  only  on  a  yearly  basis.   The  system  reserve 
capacity  must,  therefore,  deal  with  averages  over  a  yearly  interval. 
The  two  types  of  annual  fixed  costs  are  described  in  further  detail  in 
the  following  sections. 

3.7.1  Reserve  Margin  Requirement 

A  reserve  capacity  is  necessary  to  ensure  that  adequate  power  is 
available  in  the  event  of  an  unscheduled  outage  at  a  plant.   The  types 
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of  components  considered  in  this  research  contribute  to  the  overall 
plant  forced  outage  rate,  denoted  by  FOR,  and  thus  the  component  should 
be  assigned  some  of  the  cost  of  the  reserve  margin  for  the  plant.   A 
method  proposed  by  Garver  (24)  was  adopted  in  this  study. 

The  addition  of  plant  n  with  a  capacity  MW  and  FOR  to  a  utility 
system  requires  some  fraction  of  MW  to  be  placed  in  reserve  status  in 
order  to  maintain  the  same  system  reliability.  Let  RC(FOR)  be  the 
increase  in  reserve  capacity  necessary  for  the  addition  of  MW.   An 
estimate  of  RC(FOR)  is 

RC(FOR)  =  SRC  *n[(l-FOR)  +  FOR  e™^8**1]  ,         (3.65) 

where  SRC  is  the  system  risk  characteristic  of  the  utility's  generating 
system.  The  system  risk  characteristic  is  a  measure  of  the  change  in 
loss-of-load  probability,  denoted  by  LOLP,  for  a  given  change  in  reserve 
capacity  and  is  a  measure  of  the  effect  on  LOLP  by  a  change  in  reserve 
capacity.   The  SRC  is  approximated  by  plotting  the  natural  logarithm  of 
the  LOLP  for  the  system  against  the  system  reserve  capacity.   A  straight 
line  is  drawn  between  the  present  system  LOLP  on  this  curve  and  a  LOLP 
which  is  e  times  larger  on  the  same  curve,  where  e  is  2.718...  The 
negative  of  the  slope  of  the  line  drawn  is  the  inverse  of  the  SRC.   A 
typical  value  for  SRC  is  between  0.2  to  0.4  of  the  unit  generating 
capacity  (12).  An  alternate  expression  used  to  estimate  SRC  is  (24) 

N 
SRC  =  I   MW(i)  FOR(i)  ,  (3.66) 

1-1 

where  N  is  the  number  of  plants  in  the  system.   The  above  expression  for 
SRC  is  only  a  "rough  approximation"  used  to  describe  the  SRC. 
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Each  component  in  a  plant  contributes  to  the  overall  FOR.   A 
portion  of  the  reserve  capacity  can  therefore  be  assigned  to  a 
particular  component  by 

RC(FORc)  s  RC(FOR+FORc)  -  RC(FOR)  ,  (3.67) 

where  F0Rc  is  the  forced  outage  rate  attributed  to  the  component.   A 
method  of  determining  FOR  must  be  found  to  determine  the  reserve 
capacity  associated  with  the  component  of  interest. 

The  standard  definitions  of  forced  outage  rate  and  plant 
availability  factors  are 


PH  -  SOH  (3.68) 

and 

.        PH  -  SOH  -  FOH 
avail        PH         '  .  (3.69) 


where  PH  is  the  number  of  hours  in  one  year  (8760  hrs) ,  SOH  is  the 
number  of  scheduled  outage  hours  in  one  year,  and  FOH  is  the  number  of 
forced  outage  hours  in  one  year.   Manipulating  the  above  equation  gives 


PH 
PH  -  SOH  =  ,  °   f    ,  .  (3   701 

1-FOR  avail  u"  /u; 


The  amount  of  time  the  plant  will  be  shutdown  during  year  i  due  to 
failure  in  a  component  of  interest  is 


I  expected  number   average  plant 
F0Hc(i)  =  favali(1)  jof  failures     i   ^downtime  per 

in  year  i         failure  in  year  i 


(3.71) 
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The  fractional  downtime  is  the  ratio  of  actual  downtime  to  the  maximum 
downtime,  T    .   Thus,  the  average  plant  downtime  can  be  written  in 
terms  of  the  average  fractional  downtime  as 


MW 
TD(i)  =  fDI(i)  Trepl  ♦  [l-fDT(i)]  T^  ^  .       (3.72) 


In  general,  the  average  fractional  downtime  is  a  weighted  average  such 

that 

i+1  ri+1 


i  df  8(f )  fDT(f )      Jt    df  g(f)  fDI(f) 

DTk    ;  ri+1  fi+l  ■      I**"' 


^  dt'    g(t')  df    g(t') 

i  Ji 

where  X  dt '  is  the  probability  of  failure  in  dt1  about  t'  and  g(t')  is 
the  probability  the  plant  is  operating  at  t ' . 

The  expected  number  of  failures  in  year  i  is  a  very  difficult 
quantity  to  estimate.   The  number  of  failures  depends  on  the  length  of 
the  downtime  for  each  failure  because  during  this  time  a  failure  cannot 
occur.   For  a  Poisson  process  with  instantaneous  replacement,  the 
expected  number  of  failures  in  time  t  would  be  Xt.   As  an  upper  bound 
for  the  FOH  estimation,  the  value  of  {X}  failures  is  used  to 
approximate  the  expected  number  of  failures  in  any  year. 

The  forced  outage  hours  for  the  component  can  now  be  written  as 


?0}iCW    ■  favail(1)  U>  WIW"  +  ^-W*)!  t)  •  <3-74> 


Combining  the  above  with  Eqs.  (3.65),  (3.66)  and  (3.67)  gives  the 
forced  outage  rate  due  to  the  component  in  year  i  as 
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FOH    (i)  ,  MW,        i 

F0Rc(i)   "  PirhoiT  "   {X}Trepl[?DT(1)    +   ^V1"  tFJ  »"«*tl)l.    (3.75) 


The  reserve  margin  cost  due  to  the  component,  denoted  by  RC,  can  now  be 
estimated  by 


RC[FOR  (i)]  =  RC[FOR  (i)  +  FOR(i) ] 

f  FOR   (i)(eMW/SRC-l)    , 

=  SRC   In   1   ♦  £ ^j-^ (3.76) 

1  l+FOR(i)(eMW/i>      -1)J    • 


For  Plan  A,  the  average  downtime  is  the  fixed  repair/replacement 
time,  I   -,   The  forced  outage  rate  due  to  the  component  for  Plan  A 
becomes 

MW  {X}  T    [l-FOR(i)] 
F0Rc(1)  "  MTPI "  {X   Trepl}  tl-FOR(i)]     (3.77) 


making  the  reserve  margin  cost  due  to  the  component 

,  r       {XT    }[l-FOR(i)](eMW/SRC-l)1 

RC  (F0R*(1))  -  SRC  In  1  +  2£i ^ .   (3.78) 

A    °  [  l+FOR(i)(eMW/SRC-l)     I 


For  Plan  B,  the  average  fractional  downtime  is  found  by  numerically 
evaluating  Eq.  (3.75).   For  the  case  of  a  constant  fractional  downtime, 
and  g(t)  the  forced  outage  rate  due  to  the  component  is 

a  r  MW,   i 

FOR*(i)  =  {X  Trepl}  [fDT(i)  +  [l-fDT(i)]  -SJ2S.J  [l-FOR(i)]    (3.79) 

making  the  reserve  margin  cost  due  to  the  component 

MW 


RC  (FOR  (i))  =  SRC  in 


ftWlV')t'1V"-rl'1-"""w,w-" 


[l+FOR(l)(.'W(n)/S,,C-1)] 


.(3.80) 
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The  reserve  margin  cost  is  an  expenditure  that  must  be  converted  to 
a  revenue  requirement  using  the  fixed  charge  rate.   Therefore,  the 
reserve  margin  cost  associated  with  the  component  at  plant  n  in  year  m 
is 

C   (n,m)  =  RC(n,m)  FCR(n)  £   (n,m)  ,  (3.81) 

nn  rm 

where  C^fn.m)  is  the  cost  per  unit  reserve  margin  for  plant  n  in  year 
m.  An  allowance  is  made  for  escalation  in  reserve  margin  per  unit  cost 
by 


Crm(n'm)  *  Wn'0)   V1+Irm(J)!   ■  <3-82> 

where  I (m)  is  the  escalation  rate  in  year  m. 

rm  J 

Because  the  reserve  margin  is  contracted  yearly,  the  plant  must  pay 
for  the  entire  year.  The  annual  reserve  margin  costs  are  therefore 
discounted  for  the  entire  year  in  which  the  total  estimation  is  being 
performed. 

3.7.2  Spare  Purchase  Cost 

The  spare  purchase  cost  includes  the  capital  cost  of  the  component, 
transportation,  initial  storage,  and  initial  administrative  costs.  The 
spare  purchase  cost  is  shared  by  all  members  of  the  pool  and  is  a 
capitalized  cost.  The  purchase  cost  for  plant  n  is  given  by 

C    (n,0)  =  F(n,0)  5   (0)  REV(n)  ,  (3.83) 

spr  spr 

where  C   (0)  is  the  total  initial  purchase  cost  of  the  spare  to  be 
capitalized. 
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3.7.3  Spare  Storage  Cost 

An  additional  cost  associated  with  maintaining  the  spare  component 
and  managing  the  spare  pool  is  a  Plan  B  cost  only.   In  general  this  cost 
will  be  very  small  but  is  included  for  completeness.  A  fraction  of  this 
cost  is  allocated  to  each  member  of  the  pool  based  on  the  plant 
allocation  factor.   The  spare  storage  cost  in  year  m  for  plant  n  is 


Cgg(n,m)  =  F(n,m)  5ss(n,0)    H  [1+1   (j)]|  ,  (3.84) 


where  C   (n,0)  is  the  spare  storage  cost  in  year  0  for  plant  n  and 
I   (m)  is  the  escalation  rate  in  year  m. 

3.7.4  Spare  Salvage  Value 

At  the  end  of  pool  lifetime  it  is  possible  for  the  spare  component 
to  have  some  salvage  value.   The  salvage  value  at  the  end  of  pool 
lifetime  is  capitalized  by  multiplying  by  the  revenue-to-expense 
conversion  factor.   The  spare  salvage  value  is  shared  among  the  members 
in  the  pool  and  is  a  Plan  B  cost  given  by 

Wft«W  •  ^'"max5  W^^  REV(n)  '        (3-85) 

where  C_„_(n,M )  is  the  capitalized  salvage  value  of  the  spare  at  the 

end  of  pool  lifetime,  denoted  by  M   .An  allowance  is  made  for  an 

max 

escalation  in  cost  of  the  spare  salvage  value  by 


?ssv(n'm)  ■  ^sv(n'0)   "  [1+ISsv(j)]l  • 


where  I    (m)  is  the  escalation  rate  in  year 
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3.7.5  Present  Worth  of  Annual  Fixed  Costs 

Because  the  annual  fixed  costs  are  present  worthed  at  the  end  of  a 
year.  The  present  worth  of  the  annual  fixed  costs  is  estimated  by 
summing  over  the  interval  (0,1c)  where  k  is  an  integer  value,  such  that 
k  -  1  <  t  <  k.   In  general,  the  value  of  k  will  equal  the  pool  lifetime. 
The  present  worth  of  the  annual  fixed  costs  for  Plan  A  is  obtained  by 
summing  all  costs  over  the  interval  (0,t),  as  given  by 


k    C   (n,m) 

PWAFC.(n.t)  =  J     — 

A        m=l  ' 


3I1ti+*d0J>] 


I       RC  (n.m)  FCR(n)  C   (n,0)  I      I   (*)]  ,        (3.87) 
m=l  l  A  rm      £=l   m   J 


where  I   (m) 

J-l  [      1+id(J) 


,1+Irm(J>  , 

I — 7—r         -I  and  k  =  1,2,3  ...  such  that  k-1  <  t  <  k. 


The  present  worth  of  the  Plan  B  annual  fixed  costs  over  the 
interval  (0,t)  and  a  pool  lifetime  of  M    is  given  by 


PWAFCB(n,t)    =     I 
m=l 


C      (n,m)    +  C      (n,m) 
rm ss 


jljUH^W] 


C        (n,M       ) 
ssv         max 


+  C  (n,0) 

spr 


j5iti+id(J)] 


I  |RC    (n.t)    FCR(n)    C      (n,0)      I      I      (i) 


A=l 


+  F(n,0)    C      (n,0)      I      I      (t) 
SS  £=1      ss 


F(n,0)    ?        (n,0)    REV(n)    I        (M        ) 
ssv  ssv     max 


+  F(n,0)  C   (0)  REV(n)  , 
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(3.88) 


where 


I   (m)  =  B 
J-l 


1+Iss(J) 


l+id(j) 


M 
max 

i   (m  )  =  n 

ssv   max 


1  +  1    (j) 

SSV  J 


l+ld(j> 


and  k  =  1,2,3  ...  such  that  k-1  <  t  <  k. 

The  equations  derived  in  Section  3.4,  3.5,  3.6,  and  3.7  for  the 
Plan  A  and  Plan  B  present  worth  costs  are  used  in  the  computer  code 
KSUSPARE  (see  Appendix  B) .   The  total  present  worth  costs  for  Plan  A  and 
Plan  B  are  estimated  by  Eq.  (3.16)  given  in  Section  3.3.   The  estimated 
saving  is  found  by  taking  the  difference  between  the  two  plan  costs. 

3.8  Expected  Saving  and  Benefit-Cost  Ratio 

The  costs  under  the  two  managerial  plans  can  now  be  compared.   The 
present  worth  of  all  expected  costs  at  time  t  and  at  plant  n  under  the 
two  plans  is  given  by  Eq.  (3.18).   If  the  plant  is  presently  operating 
under  Plan  A  then  the  expected  saving  of  switching  to  Plan  B  for  the 
time  t  is  given  by 


SAV(n.t)  =  PWTOTAL  (n,t)  -  PWT0TAL, (n , t ) 


(3.91) 


To  obtain  the  total  expected  saving  for  all  the  utility's  plants  in  the 
pool,  the  utility  manager  would  sum  the  expected  saving  for  all  his 
plants  in  the  spare-component  pool. 

An  industry  such  as  the  electrical  utility  industry,  which  is 
regulated  by  a  government  agency  or  commission,  is  usually  required  in 
an  economic  analysis  to  calculate  an  estimate  of  the  benefits  of  a 
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certain  project  compared  to  the  cost  of  the  project  (16).  To  be 
acceptable  economically,  the  economic  benefits  must  exceed  the  cost  of 
providing  the  benefits.   Such  a  comparison  is  performed  in  a 
benefit-cost  ratio,  BCR.   The  expected  economic  benefits  of  the  spare 
pool  for  plant  n  is  the  difference  in  total  present  worth  costs  less  the 
purchase  cost  of  the  spare.  The  economic  cost  of  the  spare  pool  is  the 
initial  purchase  cost  of  the  spare.  Thus  the  BCR  for  plant  n  for  the 
time  t  can  be  expressed  as 


PWTOTAL  (n,t)  -  [PWT0TAL„(n, t)  -  C    (n,0)] 

BCR(n.t)  -  * c   (n,0)  ^ «.90) 

spr 


BCR(n,t)  =  l+f-<^,  (3.91) 

spr 


where  C    is  the  initial  capital  cost  of  the  spare  for  plant  n.   If  the 
expected  saving  is  negative,  the  benefit-cost  ratio  is  less  than  1  and 
Plan  B  should  be  rejected.   In  contrast,  if  the  expected  saving  is 
positive,  the  benefit-cost  ratio  is  greater  than  1  and  Plan  B  should  be 
accepted. 
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CHAPTER  4 
EXAMPLE  RESULTS  AND  COMPARISON  TO  PC-SPARE 

4.1  Comparison  of  KSUSPARE  With  PC-SPARE 

In  this  section  the  difference  between  the  methodology  used  in  the 
KSUSPARE  computer  code  and  the  PC-SPARE  computer  code  is  examined.   The 
reasons  for  differences  between  the  two  code  results  can  be  classified 
into  two  different  types. 

The  type  I  difference  refers  to  the  difference  in  the  treatment  of 
how  model  costs  are  incurred  and  how  they  are  brought  to  a  present 
value.   Total  present  worth  costs  in  the  KSUSPARE  code  are  determined  by 
integrating  differential  present  worth  costs  over  an  interval  of  time. 
In  contrast,  the  total  present  worth  costs  in  the  PC-SPARE  code  are 
determined  by  summing  annual  present  worth  costs  for  each  year  over  the 
time  interval.   The  KSUSPARE  code  allows  payment  of  all  failure  related 
costs  at  any  time  after  the  failure  occurs,  while  the  PC-SPARE  code 
requires  payment  of  all  failure  related  costs  along  with  yearly  costs  at 
the  end  of  each  year. 

The  type  II  difference  refers  to  the  difference  due  to  the 
fundamental  model  differences  between  the  KSUSPARE  methodology  versus 
the  PC-SPARE  methodology.  Three  major  reasons  account  for  the  type  II 
differences.   The  KSUSPARE  methodology  includes  the  probability  that  a 
plant  be  operational  before  a  failure  can  occur  in  a  differential  time 
interval  while  the  PC-SPARE  methodology  assumes  all  plants  are 
operational  at  the  beginning  of  each  yearly  interval.   The  KSUSPARE 
methodology  incorporates  the  fractional  downtime  to  estimate  how  long  the 
spare  was  used  during  the  repair/replacement  time  interval  and  the 
reserve  margin  requirement  while  the  PC-SPARE  methodology  incorporates 
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the  product  of  the  probability  of  spare  availability  and  the 

repair/replacement  time  for  the  same  estimates.   Although  not  evaluated 

in  this  section,  another  factor  contributing  to  a  type  II  difference  is 

the  use  of  a  variable  spare  availability  g  (t) .   The  exponential 

behavior  of  gg(t)  over  the  first  repair/replacement  time  interval 

increases  the  probability  a  plant  is  operating  compared  to  values  at 

later  times.   Thus,  the  likelihood  of  a  failure  in  the  first 

repair/replacement  time  increases  as  does  the  fractional  downtime.  The 

end  result  is  an  increase  in  both  plans'  total  costs.  The  difference 

caused  by  a  variable  gg  (t)  is  significant  only  when  evaluating  the 

spare-pool  costs  for  short  pool  lifetimes  (e.g.,  less  than  five 

repair/replacement  time  intervals) . 

The  following  sections  describe  how  the  type  I  difference  is 

estimated  along  with  the  type  II  difference  estimates  for:   1)  operation 

and  maintenance  costs,  2)  reserve  capacity  costs,  3)  replacement  energy 

costs,  and  4)  used  component  salvage  value  (see  Table  4.1).   To  analyze 

the  differences  between  the  computer  codes  certain  simplifying 

assumptions  were  made:  (i)  constant  component  availability  (i.e., 

constant  g(t)),  (ii)  constant  fractional  downtime,  and  (iii)  constant 

plant  capacity  factors.   For  these  assumptions  the  spare  component 

purchase  cost,  spare  final  salvage  value,  and  the  spare  storage  costs 

are  identical  for  both  studies  and  are,  therefore,  neglected  in  this 

analysis. 

Table  4.1  A  grouping  of  cost  categories  which  constitute  the  type  I 
and  type  II  differences. 


Cost  Category Type  I  difference     Type  II  difference 

Operation  and  maintenance  X  X 

Replacement  energy  X  X 

Reserve  capacity  X 

Used  component  salvage  value  X 
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The  concept  of  equivalence  factors  is  introduced  to  quantify  the 
type  I  and  type  II  differences.   The  equivalence  factor  is  defined  such 
that  if  the  i-th  cost  estimated  by  KSUSPARE  is  multiplied  by  the  i-th 
type  I  equivalence  factor,  denoted  by  X  ,  and  the  product  is  then 
subtracted  from  the  i-th  cost  estimated  by  PC-SPARE,  the  result  would 
equal  the  type  II  difference.  Thus,  the  type  I  equivalence  factor  for  a 
particular  cost  i  is  the  ratio  of  the  PC-SPARE  estimate  of  the  cost  to 
that  which  would  be  predicted  by  KSUSPARE  with  type  II  differences 
removed  (i.e.,  using  the  same  model  except  for  the  continuous  time 
involved  in  taking  the  present  worth).   The  type  II  equivalence  factor, 
denoted  by  ?,,  is  the  ratio  of  cost  i  estimated  by  PC-SPARE  to  that 
predicted  by  KSUSPARE  using  the  same  economic  method  for  obtaining  the 
present  worth  (i.e.  with  type  I  differences  removed).   Expressions  for 
these  equivalence  factors  are  derived  in  the  sections  below. 

With  these  equivalence  factors  one  can  now  quantify  the  numerical 
differences  between  PC-SPARE  and  KSUSPARE.   To  obtain  the  i-th  cost 
difference  for  plant  n,  the  term,  PWi(n, K, PC-SPARE) ,  is  used  to  denote 
the  present  worth  of  the  i-th  cost  over  K  years  estimated  by  the 
PC-SPARE  code  and  the  term,  D(n,K,j),  represents  the  j-th  difference 
between  KSUSPARE  and  PC-SPARE  after  an  integral  number  of  years  of 
spare-pool  operation.   Thus,  the  difference  between  KSUSPARE  and 
PC-SPARE  is  given  by 

D(n,K,I)  =  PWi(n,K, PC-SPARE)  -  Y±(n)  PWi(n,K, KSUSPARE)      (4.1) 
and 

D(n,K,II)  =  PWi(n,K, PC-SPARE)  -  X  (n)  PWi (n,K, KSUSPARE) .     (4.2) 

Later  in  Section  4.4  an  example  is  given  and  explicit  values  for 
differences  are  shown  for  certain  costs. 
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4.2  Equivalence  Factors  for  Type  I  Differences 

The  present  worth  of  the  1-th  cost  estimated  by  KSUSPARE  for  Plan  A 
over  the  interval  (0,K)  is  expressed  as 


PWiA(n,K,KSUSPARE)  =  qCn.O)  g   dt '  I^t'), 


(4.3) 


where 


lAf) 


m-l  fl+I^j) 

n 


jm  ira^ 


l+I1(m) 


l+id(m) 


and  m  =  INT(t')  such  that  ¥  =  t'  -  m.   Because  the  escalation  and 

m 

discount  rates  are  constant,  Eq.  (4.3)  can  be  analytically  evaluated  as 

fK  *'  aK   1 

PWOM    (n,K, KSUSPARE)    =   C      (n,0)g       dt'    aC      =   {      (n,0)    g  ~r^k,         (4.4) 
ft  om  j  _  om  Jin  (a; 


where  a  =    (1+1 . ) /(1+i.) . 
l  d 

By  contrast,  the  PC-SPARE  present  worth  of  the  i-th  cost  estimated 
by  PC-SPARE  for  Plan  A  over  the  interval  (0,K)  is  expressed  as 


PWiA(n,K, PC-SPARE)  =  C±(n,0) (1-PjI^  j)   £  I.(j) 


(4.5) 


where  Pf  is  the  probability  of  a  single  failure  in  a  one  year  interval 
equal  to  Xe   .   Equation  (4.5)  can  be  rewritten  for  constant  escalation 
and  discount  rates  as 


i        M 

PWi  (n,K, PC-SPARE)  =  % . (n,0) (l-\e~  T   ,)  J     a? 
A  i  repl   .«•, 

j-1 


=  ;,(n,0)(l-Xe"XT   ,)   1_a 


repl   K,,   . 
a  (1-a) 


(4.6) 


Thus,  the  type  I  difference  is  given  by 
:'  I,.(f) 
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df 
I  n 


0     l 

X.(n)  =  -+ ,  (4.7) 


I  I±(j) 
j=l 

where  for  constant  escalation  and  discount  rates 

_.   (a-1)  a  (1-a)  _  (a-1)  a  ,.  ., 

ln(a)  (1-a)     In (a)   "  C4-8) 

4.3  Equivalence  Factors  for  Type  II  Differences 

Type  II  differences  between  PC-SPARE  and  RSUSPARE  result  because  of 
the  different  formulations  of  the  individual  costs  in  the  two  computer 
modeling  techniques.   Thus,  the  type  II  equivalence  factor  must  be 
evaluated  separately  for  each  individual  cost.   Similar  to  the  type  I 
equivalence  factors,  the  type  II  equivalence  factors  will  be  evaluated 
over  the  time  interval  (0,K).  The  type  II  differences  for  Plan  A  and 
Plan  B  are  evaluated  for  the  1)  operation  and  maintenance  costs,  2) 
reserve  capacity  costs,  3)  replacement  energy  costs,  and  4)  used 
component  salvage  value. 

4.3.1  Operation  and  Maintenance  Costs 

The  expression  for  Plan  A  operation  and  maintenance  costs  used  by 
PC-SPARE  is 


PWOM  (n,K, PC-SPARE)  =  C   (n,0) [ l-*exp(-A  T   .)]    )   I   (j )  .     (4.9) 
A  om  r     repl   ._,   om  J 


The  expression  for  Plan  A  operation  and  maintenance  costs  used  by 
KSUSPARE  is 
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PWOMA(n,K,KSUSPARE)  =  ?   (n,0)  g   dt'  I   (t*)  .  (4.10) 

A  om       J ~  om 


Thus,  the  Plan  A  type  II  equivalence  factor  for  operation  and 
maintenance  at  plant  n  is  given  by 


1-Xexp(-X  T    ) 
Yom(n)  ■  g   "^  ■  (1+X  Trepl)[l-XeXp(-X  Trepl)].   (4.11) 


The  expression  for  Plan  B  operation  and  maintenance  costs  used  by 
PC-SPARE  is 


PW0M_(n,K, PC-SPARE)  =  C   (n,0) [ 1-Xexp(-X  T 

o  om  repl 


*exp(-XT    )]  J  I      (il    .  (4.12) 


The  expression  for  Plan  B  operation  and  maintenance  costs  used  by 

KSUSPARE  is 

K 

PW0M^(n,K, KSUSPARE)    =    C      (n,0)    g(l+XT        ,-f„   XT        .)        dt'    I      (t').      (4.13) 
a  om  repl     DT        repl      j  om 

Thus,  the  Plan  B  type  II  equivalence  factor  for  operation  and  maintenance 
is  given  by 

„  [l-Xexp(X  T        .)    +   Xexp(-X  T  )] 

YB    (n)    > - ^ - temp 

omW  g(l  +  X  T        ,-f^X  T        ,. 

repl     DT        repl) 


(1+X  T        .)[1-Xexp(-X  T        .)    +   Xexp(-X  T  )] 

repl    r  repl K  temp   ' 


(1+X  T        ,    -   f      X  T  ) 

v  repl  DT        repl/ 

4.3.2  Reserve  Capacity  Costs 

The  expression  describing  the  forced  outage  rate  due  to  the 
component  for  Plan  A  used  by  PC-SPARE  is 


(4.14) 
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FOR*(n, PC-SPARE)  =  Xexp(-X  T   . ) (l-FOR(n) ) , 
o  repl 


(4.15) 


where  FOR(n)  is  the  forced  outage  rate  of  the  plant.   In  contrast,  the 
relationship  describing  the  forced  outage  rate  due  to  the  component  for 
Plan  A  used  by  KSUSPARE  is 


FOR"(n,KSUSPARE)  -  XT   ,(l-FOR(n)) 
c  repi 


(A. 16) 


The  expression  describing  the  reserve  margin  costs  used  by  both  codes  is 


PWRM(n.K)    =  SRC(n)«.nU   + 


T0Kc(n)[emM/SRCM    -    1] 


C  (n.O)FCR(n)  •  J  I  (j)  , 
rm  .**.      rm 


(4.17) 


j   rl+I   (i) 
where     I   (j)  =  n   .  " 


Thus,  the  Plan  A  type  II  equivalence  factor  for  reserve  capacity  is 
given  by 


£n 


rm 


Xexp  (-X  Tfepl)  [  l-FOR(n)  ]  [e""M  /SRC  (n)  ] 


l+FORte™^'5^-!] 


ta 


XT        Al-F0Hn)][e™M/SRCM-l] 

1   +        "P1 


(4.18) 


1   ♦  ront.MWW/SRCCn)    _    " 


The  expression  describing  the  forced  outage  rate  due  to  the 
component  for  Plan  B  used  by  PC-SPARE  is 

MW 


B  -X 

FOR  (n.PCSPARE)  ■  X  e  (T    -P      T    +P 


down 


repl  avail  temp  avail  MW(n>   temp 


T    >[l-FOR(n)]  ,    (4.19) 
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where  P    .  is  the  probability  of  spare  availability  estimated  by  the 
PC-SPARE  code  using  binomial  theory.   In  contrast,  the  relationship 
describing  the  forced  outage  rate  due  to  the  component  for  Plan  B  used 
by  KSUSPARE  is 

FORJn.KSuSPARE)  -  X  T^  fm  ♦  (1-f^)  ^^  [  l-FOR(n)  ]  ,     (4.20) 


where  f   is  the  average  fractional  downtime  and  MW.,    is  the  derated 
DT  down 

capacity  for  plant  n  during  the  time  the  spare  is  installed. 

Thus,  the  Plan  B  type  II  equivalence  factor  for  reserve  capacity 
at  plant  n  is  given  by 


X  e  [I   -P     I   +P 


Jin  -U 


repl  avail  temp  avail  MW(n)   temp 


dOWn  T    l[l-F0R(n)I(eMW(n>/SRC(n)-  1) 


1  .  FOR(n)  (.■WWW  -13 


Y  (n) 


£n 


A  T   ,  A      ♦  Cl-i  >  ""down]  [l-F0R(„,]teMW(n)/SRC(n)  -11 
repl   M      DT  JKnT 


1  ♦  F0R(„)  [.""WSRCCn)  .  y 


.(4.21) 


4.3.3  Replacement  Energy  Costs 

The  expression  describing  the  replacement  energy  costs  for  Plan  A 
used  by  PC-SPARE  is 


PWENA(n,K, PC-SPARE)  =  Xexp(-XT    )f    (n)MW(n) ?   (n,0)   £  I   (j),   (4.22) 


repl  cap 


i-i 


wh 


j   f1+Ien(i) 


The  expression  used  by  KSUSPARE  is 


PWEN  (n.K.KSUSPARE)  =  XT   .  f    (n)MW(n) 5   (n,0)  g   dt '  I   (t').   (4.23) 
A  repl  cap        en       in  en 


in 


Thus,  the  Plan  A  type  II  equivalence  factor  for  replacement  energy  at 
plant  n  Is  given  by 

-X 

1 

"en~"'    g    N*     "repl 


*i»  -  V  =  (1  +  x  T   _)  e"X  .  (4.24) 


The  expression  describing  the  replacement  energy  for  Plan  B  used  by 
PC-SPARE  is 

PWEN  (n,K, PC-SPARE)  =  X  e~A  f    (n)  MW(n)  C   (n.ofl   . 

B                    cap  en     {   repl 

MW  i   K  . 

'  Pavail  Ttemp  +  Pavail  MW(n)   TtempJ   |  r«n^  '  (4-25) 

The  expression  used  by  KSUSPARE  is 


PWEN  (n.K.KSUSPARE)  =  X  T   .  f    (n)  5   (n,0)  g  [MW,    +  MW   f   ] 
B  repl  cap     en      •  '   down     spr  DTJ 

K 

df  I   (f)  ,  (4.26) 

0 


where  MW    is  the  spare  capacity. 

Thus,  the  Plan  B  type  II  equivalence  factor  for  replacement  energy  is 
given  by 

-X  mj 

-„       (1+AT   ,)eAMW(n)[T    -P   „    T    +P        down 

vB  fnl  .      repl repl  avail   temp  avail  MW(n)    tempJ     ,,  „„ 

en  —  •  (4.27) 

T   ,  [MW,    +  MW    (-_] 
repl   down     spr  DT 

4.3.4  Used  Component  Salvage  Value 

The  expression  describing  the  used  component  salvage  value  used  by 

PC-SPARE  is 

K  « 
PWUSV(n,K, PC-SPARE)  =  5   (n)  REV(n)   J   I    (1)  ,  (4  28) 

usv  .L,       usv  J  n,'°J 

j-1 
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m      rl+I        (i) 

where  I        (j)    =     I      L      "ff. 
usv  -"        ,_,     [1+i    (i) 


i=l  *■   d 


The  expression  for  the  used  component  salvage  value  used  by  KSUSPARE  Is 

K  . 
PWUSV(n,K, KSUSPARE)  =  X.        (n,0)  REV(n)  g  \      I    (j)  .      (4.29) 


Thus,  the  equivalence  factor  for  used  component  salvage  value  at  plant  n 
is  given  by 

Y   (n)  =  -  =  (1  +  X  T   ,)  .  (4.30) 

usv      g  repl  v.-1"/ 

4.4  Example  Problem  Comparing  KSUSPARE  and  PC-SPARE 

An  example  problem  was  designed  to  compare  the  KSUSPARE  code 
results  to  the  PC-SPARE  code  results.   The  Plan  A  and  Plan  B  costs, 
expected  savings,  and  benefit-cost  ratios  were  estimated  for  one  plant 
and  a  15  year  pool  lifetime.  Table  4.2  provides  a  listing  of  the 
problem  input  data  with  all  input  dollar  values  listed  as  present  worth 
values.   The  fractional  downtimes  estimated  by  the  three  simulation 
failure  models  as  well  as  renewal  theory  were  used.   In  addition,  the 
asymptotic  and  time-dependent  g(t)  were  used  to  estimate  the  total 
plans'  costs.   Table  4.3  includes  a  listing  of  the  computer  output  for 
both  KSUSPARE  and  PC-SPARE.   A  close  examination  of  this  table  indicates 
that  the  Plan  A  costs  estimated  by  KSUSPARE  are  larger  than  the  Plan  A 
costs  estimated  by  PC-SPARE.   The  Plan  B  costs  estimated  by  KSUSPARE, 
except  for  the  renewal  theory  estimates,  are  smaller  than  the  Plan  B 
costs  estimated  by  PC-SPARE. 
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Table  4.2  Input  data  for  sample  calculation. 


Plant  data 


discount  rate  13 

fixed  charge  rate  ,174 

time  period  for  fixed  charge .rate  (yr)  25 

cost  of  reserve  capacity  (k$  /MW)  186 

system  risk  characteristic  (MW)  1000 

Plant  component  data 

component  capacity  (MW)  800 

cost  of  replacement  component  (k$)  2770 

differential  cost  for  replacement  energy  (k$/MWY  )  74.898 

shipping  and  installation  cost  for  using  spare  (k$)  40 

annual  operation  and  maintenance  cost  (k$)  200 

component  failure  rate  (yr  )  ,01 

salvage  value  of  used  component  (k$)  400 

salvage  value  of  failed  component  (k$)  40 
repair/replacement  time  for  failed  component  (yr)  .7397(270d) 
maximum  time  spare  can  be  used  as  temporary  substitute  (yr)   .6164(225d) 

Spare  component  data 

spare  capacity  (MW)  500 

initial  cost  of  spare  component  (k$)  2200 

final  salvage  value  at  end  of  pool  life  (k$)  1000 

annual  maintenance  and  storage  (k$)  20 

spare  failure  rate  (yr~  )  _0i 

spare  pool  lifetime  (yr)  15 

Annual  plant  rates 

plant  capacity  factor  0.390 

plant  forced  outage  rate  0.220 


Annual  escalation  rates 


repaired/replacement  component      0.050     0.070     0.070 
replacement  energy  0.040     0.040     0.040 


reserve  margin  0.050  0.050  0.050 

operation  and  maintenance  0.030  0.030  0.030 

spare  storage  0.030  0.030  0.030 

spare  shipping  and  installation  0.035  0.035  0.035 

salvage  values  0.055  0.055  0.055 


Years 

0-5 

6-10 

0.050 

0.070 

0.040 

0.040 

0.050 

0.050 

0.030 

0.030 

0.030 

0.030 

0.035 

0.035 

0.055 

0.055 

11  -  15 


k$  =  $1000 

MWY  =  mega-watt  year 
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4. 4.1  Type  I  and  type  II  Differences 

One  reason  for  the  difference  in  cost  estimates  can  be  explained  by 
the  type  I  and  type  II  differences  derived  previously  in  Sections  4.1.1 
and  4.1.2.   The  type  I  and  type  II  differences  associated  with  the 
operation  and  maintenance  costs,  reserve  capacity  costs,  replacement 
energy  costs,  and  used  component  salvage  were  estimated  for  the  renewal 
theory  cost  estimates  (see  Table  4.4).   The  renewal  theory  cost 
estimates  were  multiplied  by  the  type  I  equivalence  factors  evaluated 
for  the  different  costs  (see  Table  4.5).  The  percentage  difference  (now 
only  due  to  type  II  differences)  between  the  renewal  theory  cost 
estimates  and  the  PC-SPARE  cost  estimates  again  were  calculated  (see 
Table  4.4).   A  study  of  Table  4.4  indicates  the  following  results: 

1.  The  major  portion  of  Plan  A  and  Plan  B  operation  and 
maintenance  costs,  and  Plan  A  replacement  energy  costs  are  due 
to  type  I  differences . 

2.  The  Plan  B  replacement  energy  cost  difference  is  due  to  equal 
contributions  from  type  I  and  II  differences. 

3.  The  Plan  A  and  Plan  B  reserve  capacity  cost  and  used  salvage 
value  differences  are  due  to  only  type  II  differences. 

4.4.2  Fractional  downtime  estimate 

Another  reason  for  the  smaller  Plan  B  costs  using  the  simulation 
models'  estimates  of  the  fractional  downtimes  is  because  the  component 
removal  and  installation  times  are  neglected.   Therefore,  the  fractional 
downtimes  are  underestimated  by  the  simulation  models  except  when  the 
component  removal  and  installation  times  are  very  small  compared  to  the 
component  repair/replacement  times. 
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Table  4.4.   Percentage  differences  between  PC-SPARE  results  and  KSUSPARE 
results  with  the  PC-SPARE  results  used  as  the  benchmark. 
The  renewal  theory  estimate  for  fractional  downtime  was  used 
in  KSUSPARE  calculations  for  the  case  of  asymptotic  g(t) 
and  constant  fractional  downtime. 


Cost  categories  (k$  ) 


KSUSPARE  and 

PC-SPARE 

percent  differences 


Plan  A 


Plan  B 


KSUSPARE  multiplied  by 

type  I  equivalence 

factors  and  PC-SPARE 

percent  differences 


Plan  A 


Plan  B 


Operation  and  maintenance 
Reserve  capacity 
Replacement  energy 
Used  component  salvage  value 


k$  =  $1000 


+4.8 

+4.7 

0.0 

0.1 

+  1.0 

+5.4 

+  1.0 

+5.4 

+4.6 

+9.0 

+0.3 

+4.7 

-0.7 

-0.7 

-0.7 

-.07 

Table  4.5.   Type  I  equivalence  factors 
calculated  for  the  sample 
problem. 


Equivalence  factor 


Value 


om 
rm 

en 

X    (n) 
usv 


0.95438 
1.00000 
0.95908 
1.00000 
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4.5  Time  of  Payment 

In  Chapter  3  the  assumption  was  made  that  the  costs  incurred  due  to 
a  failure  are  paid  at  the  time  of  failure.  Another  possible  scenario  is 
that  the  costs  incurred  due  to  a  failure  are  paid  at  the  end  of  the 
repair/replacement  time.  The  effect  of  such  a  payment  schedule  was 
evaluted  using  the  input  data  from  the  previously  described  example 
problem  (see  Table  4.2).   Equations  for  the  present  worth  of  the 
operation  and  maintenance  costs,  replacement  energy  costs,  and  the  used 
component  salvage  value  were  modified  by  changing  the  value  of  t'  to  t' 
+  T. 

The  percentage  difference  between  the  "paid  at  the  time  of  failure" 
costs  and  the  "paid  at  the  end  of  the  repair/replacement  time"  costs  for 
the  renewal  theory  estimates  are  listed  in  Table  4.6.   A  study  of  this 
table  reveals  that  both  the  Plan  A  and  Plan  B  "paid  at  the  end  of  the 
repair/replacement  time"  costs  are  lower  than  the  "paid  at  the  time  of 
failure"  costs.   This  result  may  have  occurred  because  the  discount  rate 
is  larger  than  any  of  the  escalation  rates.   In  contrast,  if  the 
escalation  rates  were  larger  than  the  discount  rate,  the  "paid  at  the 
end  of  the  repair/replacement  time"  costs  would  be  higher. 
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Table  4.6.  Percentage  difference  between  the  "paid  at 
the  time  of  failure"  costs  and  the  "paid  at 
the  end  of  the  repair/replacement  time"  costs 
with  the  latter  results  used  as  the 
benchmark.  The  KSUSPARE  calculations  are  for 
the  case  of  constant  g(t)  and  fractional 
downtime. 


,  Percentage  differences 

Cost  categories  (k$  )  Plan  A       Plan  B 


Operation  and  maintenance  -6.6         -6.6 

Replacement  energy  -5.9         -5.9 

Used  component  salvage  value         -5.0         -5.0 

Jk$  =  $1000 

4.6  Optimal  Number  of  Spares 

Throughout  the  present  study  only  the  case  of  one  spare  in  the 
spare-component  pool  has  been  considered.   As  the  failure  rate  of  the 
component  increases  and/or  more  plants  join  the  spare  pool,  it  may  be 
necessary  to  purchase  and  store  more  than  one  spare.  The  question  to  be 
answered  would  be  how  many  spares  should  the  spare  pool  have  for  the 
optimal  benefit-cost  ratio  for  difference  component  failure  rates  and 
number  of  plants  in  the  spare  pool.   The  input  data  from  the  previous 
example  problem  (see  Table  4.2)  were  used  to  estimate  the  benefit-cost 
ratios  as  a  function  of  time  and  component  failure  rates  using  the 
KSUSPARE  code  (see  Figure  4.1).   A  study  of  the  figure  indicates  that 
indeed  some  optimal  number  of  spares  exists. 
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10  20 

Number  of  plants  in  the  spare  pool. 


30 


FIG.  4.1.   Benefit-cost  ratio  for  different  numbers  of  plants  in 
the  spare  pool  and  different  values  of  component  failure  rates. 
The  calculation  was  performed  using  renewal  theory  to  estimate 
the  fractional  downtime  with  the  input  data  from  the  sample 
problem  in  Chapter  4. 
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Chapter  5 
SUMMARY  AND  CONCLUSION 

In  this  chapter  a  summary  and  suggestions  of  areas  for  future 
research  are  presented. 

5.1  Failure  Models 

One  of  the  objectives  of  this  study  was  to  develop  a  more  refined 
failure  model  to  predict  the  economic  benefit  of  a  pooled  inventory 
management  system.   This  study  incorporated  some  techniques  used  in 
previous  studies  in  addition  to  developing  new  techniques  which  included 
using  computer  simulation  and  renewal  theory  to  predict  spare 
availability  and  plant  operational  probability. 

5.1.1  Renewal  Theory 

In  this  study  the  probability  of  a  plant  operating  at  a  given  time 
was  characterized  by  a  system  described  by  renewal  theory.  The  plant 
being  analyzed  was  not  considered  a  potential  candidate  for  failure  up 
to  the  time  the  plant  failed  and  the  economic  analysis  was  performed. 

The  simulation  of  random  Poisson  failures  with  a  constant  repair/ 
replacement  time  agreed  with  renewal  prediction.   The  asymptotic  renewal 
theory  result  was  also  a  good  approximation  for  the  probability  of  spare 
availability  (or  the  probability  a  plant  was  operational)  for  small 
repair/replacement  times  (i.e.,  one-half  year  or  less). 

Near  the  beginning  of  the  spare  pool  lifetime,  renewal  theory 

showed  that  the  probability  of  spare  availability  and  fractional 

downtime  varied  rapidly  and  showed  several  oscillations  before 

stabilizing  at  an  asymptotic  value  at  time  >  4  T   , .   The  maximum 

-    repl 

difference  between  the  asymptotic  renewal  equation  result  and  the 
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minimum  undershoot  value  at  T     was  highly  dependent  on  the  component 
failure  rate  (i.e.,  the  larger  component  failure  rates  had  larger 
maximum  differences).   In  contrast,  the  maximum  difference  between  the 
asymptotic  renewal  equation  result  and  the  minimum  undershoot  value  at 
Trepl  WaS  inaePendent  of  the  number  of  plants  in  the  spare  pool. 

The  transient  renewal  theory  results  are  particularly  important 
when  considering  short  pool  lifetimes  (of  the  order  of  a  few  repair/ 
replacement  times) .   For  long  pool  lifetimes  the  asymptotic  renewal 
theory  result  can  be  used  to  estimate  the  probability  that  a  plant  is 
operating.   In  both  cases  the  renewal  theory  results  are  easily 
determined  and  incorporated  in  the  economic  analysis.   Because  renewal 
theory  cannot  model  the  use  of  the  spare  for  less  than  the  entire 
repair/replacement  time,  the  fractional  downtime  calculated  using 
renewal  theory  is  the  maximum  expected  fractional  downtime. 

5.1.2  Computer  Simulation 

In  this  study  computer  simulation  was  used  to  investigate  the 
effects  of  the  spare  with  and  without  a  failure  rate  and  the  decision 
rule  for  determining  if  and  for  how  long  the  spare  was  used  by  a  plant 
during  a  failure  time.  The  decision  rule  was  considered  to  be  one  of 
two  alternatives:   (i)  the  spare  is  used  when  it  becomes  available  even 
for  very  short  periods  of  time  and  (ii)  the  spare  is  only  used  if 
immediately  available  at  the  time  of  failure  and  used  for  the  entire 
repair/replacement  interval. 

All  simulation  models  as  well  as  renewal  theory  results  for 
probability  of  spare  availability  were  identical  through  the  time 
interval  (0,  T   ^) .   This  was  expected  because  any  failure  occurring  in 
this  interval  affected  the  spare  availability.   Not  until  some  time  past 
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the  first  repair/replacement  time  did  the  spare  failure  rate  or  decision 
rule  need  to  be  considered. 

An  analysis  of  the  simulation  results  indicated  there  was  very 
little  difference  between  a  spare  with  or  without  a  failure  rate  in  a 
system  of  five  or  more  plants.   The  difference  was  negligible  in  ranges 
of  component  failure  rates  considered  by  this  study.   However,  for 
systems  with  larger  component  failure  rates  the  probability  of  spare 
availability  was  somewhat  lower  and  the  fractional  downtime  a  slightly 
larger  for  a  spare-component  pool  containing  a  spare  with  a  failure  rate 
equal  to  the  other  system  components  versus  a  spare  with  zero  failure 
rate. 

The  decision  rule  used  resulted  in  a  significant  effect  on  both  the 
spare  availability  and  fractional  downtime.   For  small  component  failure 
rates  the  decision  rule  had  little  effect  on  the  probability  of  spare 
availability,  however,  for  large  component  failure  rates  the  "use  the 
spare  only  if  immediately  available"  showed  a  much  larger  probability  of 
spare  availability.   As  expected,  the  difference  between  the  probability 
of  spare  availability  for  "use  the  spare  when  available"  and  "use  the 
spare  only  if  immediately  available"  was  much  greater  for  larger 
component  failure  rates. 

The  fractional  downtimes  were  significantly  different  for  all 
values  of  component  failure  rates.   The  combination  of  the  plants  in  the 
spare  pool  (except  the  plant  that  fails  at  the  analysis  time)  operating 
under  "use  the  spare  when  available",  a  non-zero  spare  failure  rate,  and 
the  plant  failing  at  the  analysis  time  operating  under  "use  the  spare 
only  if  immediately  available"  resulted  in  the  largest  fractional 
downtime  and  the  largest  Plan  B  costs.   The  increase  in  the  fractional 
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downtime  due  to  the  non-zero  spare  failure  rate  was  less  than  the 
decrease  in  the  fractional  downtime  due  to  the  plant  operating  under 
"use  the  spare  only  if  immediately  available".  The  fractional  downtime 
calculated  using  renewal  theory  to  estimate  the  probability  of  spare 
availability  is  therefore  considered  to  be  an  estimate  of  the  maximum 
expected  fractional  downtime. 

5.2  Economic  Analysis 

Another  objective  of  this  study  was  to  develop  an  extension  of  the 
economic  models  used  in  earlier  studies.  The  refined  probability 
failure  models  were  incorporated  into  the  generalized  economic  model  to 
obtain  the  general  cost  of  the  two  component  management  plans.  The 
economic  costs  were  also  estimated  on  a  continuous-time  basis  instead  of 
the  typical  year-end  basis. 

5.2.1  Comparing  KSUSPARE  to  PC-SPARE 

The  failure  dependent  variable  costs  (i.e.,  operation  and 
maintenance  costs)  were  a  significantly  larger  proportion  of  the  total 
plan  costs  for  KSUSPARE  while  the  difference  between  the  Plan  A  and  Plan 
B  costs  remained  approximately  the  same  as  the  PC-SPARE  results.   Using 
the  simulation  results  for  fractional  downtime  resulted  in  significantly 
smaller  Plan  B  costs  because  the  simulation  models  did  not  include  any 
allowance  for  removal  or  installation  times.   As  expected,  the 
fractional  downtime  estimated  by  renewal  theory  provided  a  maximum  Plan 
B  cost  and  the  smallest  benefit-cost  ratio. 

The  major  reasons  for  the  differences  in  costs  between  KSUSPARE 
results  and  PC-SPARE  results  were  due  to  the  method  used  to  estimate  the 
fraction  of  the  repair/replacement  time  the  failed  plant  would  be 
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shutdown  or  derating  and  the  treatment  of  how  model  costs  were  incurred 
and  how  they  were  brought  to  a  present  value.   In  addition,  a  comparison 
of  the  differences  for  the  case  of  payment  at  the  end  of  the 
repair/replacement  time  after  a  failure  versus  payment  at  the  beginning 
of  the  repair/replacement  time  showed  a  significant  decrease  in  some 
costs. 

Evaluating  the  renewal  equation  numerically  or  assuming  an 
asymptotic  value  resulted  in  no  significant  difference  in  the  estimated 
savings  or  benefit-cost  ratios  for  long  pool  lifetimes.   The  total  Plan 
A  and  Plan  B  costs  were  different  but  the  savings  between  Plan  B  and 
Plan  A  were  approximately  the  same.   For  this  reason  the  extra  expense 
and  time  necessary  to  evaluate  the  renewal  equation  numerically  is  not 
justified  unless  either  very  short  pool  lifetimes  or  very  large 
repair/replacement  times  are  being  considered. 

5.2.2  Advantages  of  KSUSPARE 

A  major  improvement  of  KSUSPARE  is  the  ability  to  handle  repair/ 
replacement  times  of  greater  than  one  year.   In  addition,  costs  incurred 
due  to  a  failure  can  be  paid  at  any  time  in  the  future  after  the  failure 
instead  of  at  the  end  of  the  year  in  which  the  failure  occurs.   Another 
improvement  is  the  ability  to  estimate  costs  as  a  continuous  function  of 
time  instead  of  at  the  end  of  each  year. 

5.3  Suggestions  for  Further  Study 

A  logical  step  would  be  to  determine  the  optimal  number  of  spares 
in  the  pool.   As  the  failure  rate  of  the  component  increases  and/or  more 
plants  join  the  pool  it  may  be  necessary  to  purchase  and  store  more  than 
one  spare  to  maintain  a  positive  savings  of  Plan  B  over  Plan  A.   Because 
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the  benefit-cost  is  a  function  of  the  spare  availability,  an  optimal 
number  of  spares  is  necessary  to  ensure  that  the  spare  is  not  always 
being  used  when  needed.   Another  area  of  study  might  be  an  attempt  to 
quantify  the  variance  associated  with  the  benefit-cost  ratio.   Of 
concern  in  such  probabilistic  studies  is  not  only  the  expected  cost  but 
the  variation  in  the  expected  cost.   The  analysis  of  such  a  variation 
would  be  difficult  to  do  with  renewal  theory  and  analytical  models  but 
easy  to  incorporate  in  simulation  models  using  Monte  Carlo  techniques. 

In  further  studies  increased  realism  could  be  added  to  the  analysis 
by  incorporating  a  yearly  load  following  curve.   The  addition  of 
seasonal  variations  in  demand  for  electricity  would  be  more  realistic 
than  the  yearly  capacity  factors  presently  used.   The  escalation  and 
discount  rates  could  also  be  described  as  continuous  functions  of  time. 
Increased  realism  could  also  be  added  by  including  the  interdependence 
of  a  failure  at  one  plant  upon  the  forced  outage  rates  and  resulting 
reserve  capacity  requirements  for  other  plants  in  the  pool. 
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Appendix  A:   SIMULATION  Computer  Code 

The  computer  code,  SIMULATION,  was  written  in  FORTRAN  and  consists  of 
a  main  program  and  10  subroutines.   A  large  number  of  variables  were 
placed  in  "common"  instead  of  being  "passed"  when  a  particular 
subroutine  was  called.   This  variable  transfer  method  resulted  in 
substantial  computer  run  time  savings. 

The  main  program  calls  a  number  of  subroutines  that  perform  the 
actual  analysis.   A  subroutine  to  read  in  the  code  parameters  is  first 
called.   The  next  subroutine  called  simulates  the  time  intervals  by 
following  the  method  outlined  in  Section  2.5  which  requires  simulating 
failures  by  calling  the  random  number  generating  subroutine  RANDOM. 
Subroutine  RANDOM  returns  a  random  number  on  the  interval  (0,1)  based  on 
an  input  seed  value.   The  seed  value  is  multiplied  by  a  constant 
followed  by  a  modulus  operation  with  a  large  number  A.   The  result  is 
then  divided  by  another  large  number  (B  which  is  A+l)  and  set  equal  to 
the  random  number.   The  next  seed  value  used  is  the  result  obtained  from 
the  modulus  operation. 

The  probability  of  failure  in  the  first  interval  (0,  TREPL)  and  the 
system  characteristic  failure  rate  are  calculated  next.  The  failure 
times  are  then  simulated  by  calling  either  TIMES  1  or  TIMES3.   After  each 
entire  simulation  is  completed  one  of  the  CHECK  subroutines  is  called  to 
determine  the  sums  used  in  subroutine  RESULT  to  estimate  the 
probability  of  spare  availability  and  fractional  downtime.   After  all 
the  simulations  have  been  run  then  the  results  are  calculated. 

Subroutine  TIMES  1  simulates  failures  using  the  process  outlined  for 
Model  1  and  Model  2.   Failures  are  simulated  until  the  last  failure 
simulated  exceeds  the  value  of  TMAX.   The  entire  simulated  time  scheme 
is  then  transferred  to  the  CHECK1  subroutine  for  analysis. 
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Subroutine  TIMES3  simulates  failures  using  the  process  outlined  for 
Model  3.  Like  in  subroutine  TIMES1,  the  entire  time  scheme  is  simulated 
and  then  the  results  are  transferred  to  subroutine  CHECK  3. 

The  total  number  of  times  a  failure  occurs  in  the  interval  (t  - 
6,  T  )  for  each  value  of  T.  is  determined  for  use  in  the  Model  1 
probability  of  spare  availability  and  fractional  downtime  estimations 
made  in  subroutine  RESULT.   In  addition,  the  total  number  of  times  and 
difference  in  time  between  T  and  the  last  time  the  spare  was  used  for 
each  value  of  t   is  determined  for  use  in  the  Model  2  probability  of 
spare  availability  and  fractional  downtime  estimations  made  in 
subroutine  RESULT. 

The  total  number  of  times  a  failure  occurs  in  the  interval  (t   - 
8i  t^)  and  difference  between  the  repair/replacement  time  and  the  actual 
downtime  is  determined  for  each  value  of  t   for  use  in  the  Model  3 
probability  of  spare  availability  and  fractional  downtime  estimations 
made  in  subroutine  RESULT.   The  last  subroutine,  PRINT3,  prints  out  the 
probability  of  spare  availability  and  its  associated  standard  deviation, 
and  the  fraction  downtime. 

Input  Data 

TREPL  =  the  failed  component  procurement  time. 

NP  =  the  number  of  plants  in  the  pool. 

LAMBDA  ■  the  component  failure  rate. 

MAXIT  =  the  number  of  simulations  the  code  is  run  for. 

FSEED  =  the  initial  seed  value  for  the  random  number 
subroutine. 

METHOD  =  determines  the  mode  which  the  code  simulates.   An 

input  value  of  1  simulates  models  1  and  2;  an  input 
value  of  2  simulates  Model  3. 
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TMAX  =  the  maximum  evaluation  time  at  which  probability  of 

spare  availability  and  fraction  downtime  are  estimated. 

NIBTR  =  the  number  of  time  intervals  between  0  and  TREPL  that 
determine  the  evaluation  times. 

NTMAX  =  the  number  of  time  intervals  between  TREPL  and  TMAX 
that  determine  the  evaluation  times. 

NREAD  =  determines  the  evaluation  times.   An  input  value  of 
other  than  999  results  in  the  evaluation  times  being 
read  in  as  input  data.   An  input  value  of  999  calls 
subroutine  GPROG. 


Code  Listing 

The  following  is  a  listing  of  the  computer  code  SIMULATION. 
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C«**«**«**«**««««»*««*******SiMULATION********************************* 

C» 

C* 

C*  THIS  PROGRAM  ESTIMATES  THE  PROBABILITY  OF  BEING  IN  PLAN  B  OR 

C«  THE  FRACTION  DOWN  TIME  BY  SIMULATING  THE  TIME  INTERVALS  BETWEEN 

C*  FAILURES  USING  NEGATIVE  EXPONENTIAL  DISTRIBUTIONS.  THE  KEY 

C»  ASSUMPTION  FOR  MODELS  1,  2,  AND,  3  IS  THE  COMPONENTS  AND  THE 

C*  SPARE  ARE  REMOVED  AND  INSTALLED  INSTANTANEOUSLY.  THE  PROB- 

C*  ABILITY  OF  BEING  IN  PLAN  B  AND  THE  FRACTION  DOWN  TIME  ARE 

C*  ESTIMATED  IN  EACH  MODEL.  ALL  TIMES  ARE  DIMENSIONLESS  (I.E., 

C*  REAL  TIME  «  COMPONENT  FAILURE  RATE). 

C« 

C*  THE  PROGRAM  USES  THREE  MODELS: 

C* 

C*  MODEL  1 

C»      

C»  THE  FAILURE  RATE  OF  THE  SPARE  IS  EQUAL  TO  ZERO.  PLANTS  IN  THE 

C*  POOL  AND  THE  FAILED  PLANT  USE  THE  SPARE  WHEN  IT  BECOMES  AVAIL- 

C*  ABLE,  EVEN  FOR  VERY  SHORT  PERIODS  OF  TIME. 

C* 

C»  MODEL  2 

C»      

C«  THE  FAILURE  RATE  OF  THE  SPARE  IS  EQUAL  TO  ZERO.  IF  THE  SPARE 

C*  IS  NOT  AVAILABLE  AT  THE  TIME  OF  A  PLANT  FAILURE,  THE  PLANT 

C*  REMAINS  SHUTDOWN  UNTIL  THE  REPAIR/REPLACEMENT  COMPONENT  IS 

C*  INSTALLED.  SIMILARLY,  IF  THE  SPARE  IS  NOT  AVAILABLE  AT  THE 

C*  TIME  OF  THE  PLANT  FAILURE  IT  WILL  NOT  BE  INSTALLED  AT  A  LATER 

C»  DATE. 
C* 

C»  MODEL  3 

C*      

C»  THE  SPARE  HAS  A  FAILURE  RATE  EQUAL  TO  THE  FAILURE  RATE  OF  THE 

C*  OTHER  SIMILAR  COMPONENTS  IN  THE  POOL.  PLANTS  IN  THE  POOL  AND 

C«  THE  FAILED  PLANT  USE  THE  SPARE  WHEN  IT  BECOMES  AVAILABLE,  EVEN 

C*  FOR  VERY  SHORT  PERIODS  OF  TIME. 

C» 

C» 

C 

C 
C 

c 
c 

c  ************** 

C  "MAIN  PROGRAM* 

C  »»**♦»****♦*** 

COMMON  T(100),SUM2(100),RATI02(100),TIME(100),TUSED(100),DOWN1(100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

COMMON/RNG/FSEED.U 

DOUBLE  PRECISION  FSEED.U 

CALL  READ 

CALL  ITER 
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CALL  PRINT 

STOP 

END 
C 
C 
C 
C 

Z******t*******t*********************t**t**************tt****tt********4 
C#  START  OF  SUBROUTINES  i 

C*t**********t*H***i****l******t**********i****t*********************t* 
C 
c 
c 
c 

c  ***************** 

C  //SUBROUTINE  READ* 

C  ***************** 

SUBROUTINE  READ 

COMMON  T(100),SUM2(100),RATI02(100),TIME(100),TUSED(100),DOWN1(100 
1 ) ,DOWN2( 100) ,S1 ( 100) .RATI0 1 ( 100) ,SUM1 ( 100) ,PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

COMMON/MAXT/TMAX,MAXTAU,NREAD,NIBTR 

COMMON/RNG/FSEED.U 

COMMON/REPAIR/DEV 

DOUBLE  PRECISION  FSEED.U 

C     

C     READ  IN  REPLACEMENT  TIME,  NUMBER  OF  PLANTS,  COMPONENT  FAILURE 

C     RATE,  TEST  VALUE,  MAXIMUM  NUMBER  OF  ITERATIONS,  AND  INITIAL  SEED. 

C     

READ    (5,10)    BETA,NP,MAXIT,FSEED 
10  FORMAT    C/,F10.6,I5,I10,D15.8) 

C  

C     READ  IN  THE  SIMULATION  MODEL  TO  BE  USED. 
C     


READ  (5,20)  MODEL 
20  FORMAT  (/,I5) 


READ  IN  THE  MAXIMUM  TIME 
FRACTION  DOWN  TIME  AT. 

TO 

EVALUATE 

THE  PROBABILITY  AND  OR 

THE 

READ 

(5, 

,20) 

MAXTAU 

READ 

IN 

THE 

NUMBER 

OF 

TIME 

INTERVALS 

BETWEEN  0 

AN! 

I  BETA. 

READ 

(5 

,20) 

NIBTR 

READ 

IN 

THE 

TOTAL  NUMBER 

OF 

TIME  VALUES  BETWEEN  0 

AND  TMAX. 

READ 

(5, 

,20) 

NTMAX 

READ 
ERROF 

IN  THE  NUMBER 
!  ESTIMATE. 

OF 

STANDARD  DEVIATIONS  USED 

TO 

CALCULATE 

THE 

READ  (5,50)  DEV 
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50  FORMAT  (/,F10.3) 
C     

C     DETERMINE  IF  THE  TIME  VALUES  TO  EVALUATE  THE  PROBABILITY  ARE  READ 
C     IN  OR  CALCULATED  WITH  A  GEOMETRIC  PROGRESSION  METHOD.  IF  THE  VALUE 
C     IS  999  THEN  THE  PROGRESSION  METHOD  IS  USED  OTHERWISE  THE  TIME 
C     VALUES  MUST  BE  READ  IN. 

C     

READ  (5,20)  NREAD 

IF  (NREAD. NE. 993)  GOTO  80 

C     

C     READ  IN  THE  TIME  VALUES  TO  EVALUATE  THE  FRACTION  DOWN  TIME  AT. 

C     

DO  30  J=1,NTMAX,1 

READ    (5,40)    TIME(J) 
40  FORMAT   (F10.4) 

S1(J)=0.0 
SUM1(J)=0.0 
SUM2(J)=0.0 
30  CONTINUE 

NTMAX2=NTMAX 
GOTO   70 
80  CALL  GPROG 
70  RETURN 
END 

C 

C 

C 

C 

C  ****************** 

C  0SUBROUTINE  GPROG* 

C  ****************** 

SUBROUTINE  GPROG 

COMMON  T( 100) ,SUM2( 100) ,RATI02( 100) ,TIME( 100) ,TUSED( 100) ,D0WN1 (100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

COMMON/MAXT/TMAX .MAXTAU , NREAD , NIBTR 

TMAX=MAXTAU»BETA 

NTMAX2=NTMAX+1 

C  

C            THE  FIRST  TIME   INTERVAL  IS  GIVEN  THE  VALUE    .01    TIMES  THE  VALUE  OF 
C            THE  FIRST  EQUALLY  SPACED   VALUE.    THE  NEXT  NIBTR-2   INTERVALS  ARE 
C            EQUALLY   SPACED.    THE  NIBTR  TIME   POINT  IS  CHOSEN   AS  EQUAL  TO  BETA. 
C  

TINT=BETA/NIBTR 

TIME(1)=.01»TINT 

NBT=NIBTR+1 

DO  20  K=2,NBT,1 

TIME(K)=(K-1)»TINT 
20  CONTINUE 

TINT2= ( TMAX-BETA ) / ( NTMAX-NIBTR ) 

C  

C  THE  EVALUATION  TIMES  BETWEEN  TMAX   AND   BETA  ARE  EITHER   EQUALLY 

C  SPACED  OR  DIVIDED  BY   A  GEOMETRIC   PROGRESSION  TECHNIQUE. 

C  
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IF  (NREAD.EQ.99D    THEN 

C  

C  TECHNIQUE  FOR  EQUALLY   SPACED  TIMES. 

C  


DO   30  K=  NBT,NTMAX,1 

TIME(K+1)=TIME(K)+TINT2 
30     CONTINUE 
ELSE 

C         

C         THE  TIME  INTERVALS  BETWEEN  BETA  AND  TMAX  ARE  DETERMINED  BY  A 
C         GEOMETRIC  PROGRESSION  METHOD  WHERE  THE  COMMON  RATIO  R  IS  THE 
C         RATIO  OF  INTERVAL  WIDTH  TO  THE  PREVIOUS  INTERVAL  WIDTH.  THE 
C         COMMON  RATIO  MUST  BE  A  VALUE  GREATER  THAN  1  FOR  THIS  METHOD. 

C         

R  =  1 . 1 

A1=(TMAX-BETA)*(1-R)/(1-R»«(NTMAX-NIBTR)) 

TIME(NIBTR+2)=A1+TIME(NIBTR+1) 

NBR=NIBTR+3 

DO    10  J=NBR,NTMAX2,1 

TIME(J)=TIME(J-1)+A1»R*»(J-NBR+1) 
10     CONTINUE 
END  IF 

C     

C     ALL  SUMS  ARE  ZEROED  IN  THIS  SUBROUTINE  FOR  THE  CALCULATED  TIMES. 

C     

DO  40  J=1,NTMAX2,1 
S1(J)=0.0 
SUM1(J)=0.0 
SUM2(J)=0.0 
40  CONTINUE 
RETURN 
END 
C 
C 
C 
C 

C  ***************** 

C  #SUBROUTINE  ITER# 

C  ***************** 

C 

SUBROUTINE  ITER 

COMMON  T( 100) ,SUM2( 100) ,RATI02( 100) ,TIME( 100) ,TUSED( 100) ,DOWN1 ( 100 
1),DOWN2(100),S1(1O0),RATIO1(1O0),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

COMMON/ERROR/RK 100) ,R2( 100) ,R3( 100) ,R4( 100) ,SA1 ( 100) ,SA2( 100) ,SA3 
1(100), 344(100) 
COMMON/RNG/FSEED,U 
COMMON/REPAIR/DEV 
INTEGER  Z,K,H,NP,MAXIT 
DOUBLE  PRECISION  FSEED,U,X 

C     

C     CALCULATE  THE  CONSTANT  TO  DETERMINE  IF  A  FAILURE  OCCURS  IN  THE 
C     INTERVAL  OF  LENGTH  EQUAL  TO  THE  REPLACEMENT  TIME  BETA  FOR  USE 
C     IN  THE  SUBROUTINE  TIMES  1. 
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PFTR=1-EXP(-(NP-1)»BETA) 
C     

C     THE  SIMULATION  STARTS  HERE. 
C     

DO  10  NIT  =  1,MAXIT,1 

IF  (MODEL. EQ. 3)  GOTO  3 
CALL  TIMES1 
CALL  CHECK  1 
GOTO  10 
3     CALL  TIMES3 
CALL  CHECK 3 
10  CONTINUE 
CALL  RESULT 
RETURN 
END 
C 
C 
C 
C 

C  #############*##### 

C  #SUBROUTINE  TIMES 1# 

C  ################### 

SUBROUTINE  TIMES 1 

COMMON  T( 100) ,SUM2( 100) ,RATI02( 100) ,TIME( 100) ,TUSED( 100) ,D0WN1 ( 100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 
COMMON/RNG/FSEED.U 

DIMENSION  DELTA* 100) ,NFAIL( 100) ,DT( 100) 
DOUBLE  PRECISION  FSEED,U,X 
INTEGER  Z,K,H,NP,MAXIT 

C     

C     DEFINE  FREQUENTLY  USED  FUNCTIONS. 
C     


FUNC1(X,J)=-DLOG(1.00000000DO-X)/J 

FUNC3(Y,X,J)=(Y-DLOG(1.00000000DO-X))/J 

FUNC4(Y,X,C)=(Y-DLOG(1.00000000DO-X))/C 

FUNC5(V,Y,C)=1.0-(EXP(-V»Y)*C) 

FUNC6(Y)=EXP(Y) 

C     

C  INITIALIZE  SUMS  AND   VECTORS  USED   IN  THE  INNER  LOOPS.    THE  VECTORS 

C  MUST  BE  REZEROED  OVER  THE  RANGE  GIVEN  IN  THE  DIMENSION  STATEMENT. 

C  

1=0 

DO  76  M=  1,100,1 
DELTA(M)=0.0 
T(M)=0.0 
DT(M)=0.0 
TUSED(M)=0.0 
76  CONTINUE 

C     

C     INCREMENT  COUNTER  TO  KEEP  TRACK  OF  THE  NUMBER  OF  TIME  INTERVALS 

C     BEFORE  TMAX  IS  EXCEEDED  FOR  EACH  SIMULATION. 

C     
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5  1=1+1 
C  

C  SET  THE  FIRST  DELTA  VALUE  EQUAL  TO  THE  REPLACEMENT  TIME,    BETA. 

C 

DELTA(1)=BETA 

C  

C     INITIALIZE  SOME  COUNTERS  AND  TIME  KEEPERS. 
C     


TUSED(I+1)=TUSED(I) 

TNEW=TUSED(I)+BETA 

F=0.0 

DSUM=0.0 

XF=0.0 

NFAIL(I)=0 

C  

C     GENERATE  RANDOM  NUMBER  BY  CALLING  ROUTINE  RANDOM. 

C     

CALL  RANDOM 

C     

C     IF  IT  IS  THE  FIRST  INTERVAL  IN  THE  SIMULATION  TRIAL  THE  FOLLOWING 
C     STEPS  ARE  BYPASSED. 

C     

110  IF  (I.EQ.1)  THEN 

DT(I)=FUNC1(U,NP) 

T(I)=DT(I) 

TUSED(I+1)=T(I) 

C         

C         CHECK  PRESENT  TIME  AGAINST  THE  MAXIMUM  EVALUATION  TIME. 

C         

IF  (T(I)-TIME(NTMAX2))  23,77,77 

C         

C         THE  LENGTH  OF  THE  INTERVAL  FOR  THE  FIRST  FAILURE  HAS  BEEN 
C         DETERMINED.  CONTINUE  IN  THE  PRESENT  SIMULATION. 

C         

23     GOTO  5 
ELSE 
END  IF 

C     

C     IF  THE  LENGTH  OF  THE  PREVIOUS  INTERVAL  IS  GREATER  THAN  BETA,  THEN 
C     (NP-1)  PLANTS  ARE  OPERATING  AT  THE  START  OF  THE  INTERVAL. 

C     

IF  (DT(I-1).GE.BETA.0R.I.EQ.2)  THEN 

C         

C  THE  VALUE  OF  THE  RANDOM  NUMBER   IS  COMPARED  TO  THE  PROBABILITY 

C  OF  FAILURE  IN  THE  INTERVAL  BETA. 

C  

IF   (U.LT.PFTR)    THEN 

DT(I)=FUNC1(U,(NP-D) 
ELSE 

DT(I)=FUNC3(BETA,U,NP) 
END   IF 

T(I)=DT(I)+T(I-1) 
IF   (T(I).GE.TNEW)    TUSEDQ  +  1 )  =T(I) 
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C         CHECK  PRESENT  TIME  AGAINST  THE  MAXIMUM  EVALUATION  TIME. 
C         

IF  (T(I)-TIME(NTMAX2))  21,77,77 
C         

C         THE  LENGTH  OF  THE  INTERVAL  FOR  THE  CASE  OF  DT(I-1)lBETA  HAS 

C         BEEN  DETERMINED.  CONTINUE  IN  THE  PRESENT  SIMULATION. 

C         

24  GOTO  5 
ELSE 

ICOUNT=0 

11=1-2 

IF  CII.GE.(NP-D)  II  =  (NP-1) 

DO  14  J=  1,11,1 

DELTA( J+1 )=BETA-T(I-1 )+T(I-J-1 ) 

ICOUNT=ICOUNT+1 

IF  (DELTAtJ+O)  13,13,11 

13  DELTA(J+1)=0.0 
ICOUNT=ICOUNT-1 
GOTO  17 

14  CONTINUE 

C         

C         THE  NEXT  LOOP  DETERMINES  WHICH  F(T)  FUNCTION  TO  USE. 

C         

17  NFAIL(I)=ICOUNT 

IF  (ICOUNT.GE.(NP-D)  THEN 
DO  11  J=1,II,1 

DELTA(J)=DELTA(J)-DELTA(II+1) 
1 1         CONTINUE 

TWAIT=DELTA(II+1) 
H=NFAIL(I)+1 
DELTA(II+1)=0.0 
ELSE 

TWAIT=0.0 
H=NFAIL(I)+2 
END   IF 

DO  75  K=  1,H,1 
ALPHA=NP-H+K 
DSUM=DSUM+DELTA(H+1-K) 
EP=FUNC6(DSUM) 
IF  (K.EQ.H)  THEN 
F=1.0 
GOTO  29 
ELSE 

F=FUNC5( ALPHA, DELTA ( H-K ) ,EP) 
END  IF 
29         XF=FUNC4(DSUM,U, ALPHA) 
IF  (U.LE.F)  THEN 
DT(I)=XF 

T(I)=T(I-1)+DT(I)+TWAIT 
IF  (T(I).GE.TNEW)  TUSED(I+1 )=T(I) 
IF  (T(I)-TIME(NTMAX2))  25,77,77 

25  GOTO  5 
ELSE 

END  IF 
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75     CONTINUE 
END  IF 
C     

C     THE  LENGTH  OF  THE  INTERVAL  HAS  NOW  BEEN  DETERMINED  FOR  THE  CASE 

C     DT(I-1)[BETA.  CONTINUE  IN  THE  PRESENT  SIMULATION. 

C     

GOTO  5 
77  RETURN 

END 
C 
C 
C 
C 

c  ******************* 

C  0SUBROUTINE  TIMES3* 

C  ******************* 

SUBROUTINE  TIMES3 

COMMON  TC 100) ,SUM2( 100) ,RATI02( 100) ,TIME( 100) ,TUSED( 100) ,D0WN1 (100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I.NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

COMMON/RNG/FSEED,U 

DIMENSION  DELTA( 100) ,NFAIL( 100) ,DT( 100) 

DOUBLE  PRECISION  FSEED,U,X 

INTEGER  Z,K,H,NP,MAXIT 

C     

C     DEFINE  FREQUENTLY  USED  FUNCTIONS. 
C     


FUNC1(X,J)=-DLOG(1.00000000DO-X)/J 
FUNC3(Y,X,J)=(Y-DLOG(1.0O0000OODO-X))/J 
FUNC4(Y,X,C)=(Y-DLOG( 1 .00000000D0-X) )/C 
FUNC5(V,Y,C)=1.0-(EXP(-V»Y)«C) 
FUNC6(Y)=EXP(Y) 

C     

C     INITIALIZE  SUMS  AND  VECTORS  USED  IN  THE  INNER  LOOPS.  THE  VECTORS 
C     MUST  BE  REZEROED  OVER  THE  RANGE  GIVEN  IN  THE  DIMENSION  STATEMENT. 

C     

1=0 

DO  76  M=  1,100,1 
DELTA(M)=0.0 
T(M)=0.0 
DT(M)=0.0 
76  CONTINUE 

C     

C     SET  LOGICAL  VARIABLE  VALUES. 
C     


IFLAG=0 
IM=1 
C     

C  INCREMENT  COUNTER  TO  KEEP  TRACK   OF  THE  NUMBER  OF  TIME  INTERVALS 

C  BEFORE  TMAX   IS  EXCEEDED  FOR   EACH  SIMULATION. 

C  

5  1=1+1 

C  

C     INITIALIZE  SOME  COUNTERS  AND  TIME  KEEPERS. 
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F=0.0 

DSUM=0.0 

XF=0.0 

NFAIL(I)=0 

C     

C     GENERATE  RANDOM  NUMBER  BY  CALLING  ROUTINE  RANDOM. 
C     

CALL  RANDOM 
C     

C     IF  IT  IS  THE  FIRST  OR  SECOND  INTERVAL  IN  THE  SIMULATION  TRIAL  THE 
C     FOLLOWING  STEPS  ARE  BYPASSED. 

C     

110  IF  (IM.LT.3)  THEN 

DT(I)=FUNC1(U,NP) 
IF  (I.EQ.1)  THEN 

T(I)=DT(I) 
ELSE 

T(I)=DT(I)+T(I-1) 
END  IF 

C         

C         CHECK  PRESENT  TIME  AGAINST  THE  MAXIMUM  EVALUATION  TIME. 

C         

IF  (T(I)-TIME(NTMAX2))  23,77,77 

C         

C         THE  LENGTH  OF  THE  INTERVAL  HAS  BEEN  DETERMINED.  THE  COUNTER 
C         IS  NOW  INCREASED  AND  THE  SIMULATION  IS  CONTINUED. 

C         

23     IM=IM+1 
GOTO  5 
ELSE 
END  IF 

C     

C     IF  THE  LENGTH  OF  THE  PREVIOUS  INTERVAL  IS  GREATER  THAN  BETA,  THEN 
C     NP  PLANTS  ARE  OPERATING  AT  THE  START  OF  THE  INTERVAL. 


C 


IF  (DT(I-I).GE.BETA)  THEN 
DT(I)=FUNC1(U,NP) 
T(I)=DT(I)+T(I-1) 
IFLAG=0 

IF  (T(I)-TIME(NTMAX2))  24,77,77 
24     GOTO  5 
ELSE 

IF  (NP.LE.1)  THEN 
1=1+1 

DT(I)=BETA+T(I-NP)-T(I-1) 
T(I)=DT(I)+T(I-1) 
IF  (T(I)-TIME(NTMAX2))  36,77,77 
36         GOTO  5 
ELSE 
END  IF 
END  IF 
ICOUNT=0 
11=1-2 
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IF  CII.GE.NP)  II=NP 

DO  14  J=  1,11,1 

DELTA(J)=BETA-T(I-1)+T(I-J-1) 

ICOUNT=ICOUNT+1 

IF  (DELTA(J))  13,13,11 

13  DELTA(J)=0.0 
ICOUNT=ICOUNT-1 
GOTO  17 

14  CONTINUE 


c 
c 

c 

17 

THE 

NEXT  LOOP  DETERMINES  WHICH  F(T)  FUNCTION  TO  USE. 

NFAIL(I)=ICOUNT 

H=NFAIL(I)+1 

IF 

(IFLAG.EQ.O)  THEN 
NN=0 

ELSE 

NN  =  1 

END 

IF 

DO 

75  K=  1,H,1 
ALPHA=NP-H+K+NN 
DSUM=DSUM+DELTA(H+1-K) 
EP=FUNC6(DSUM) 
IF  (K.EQ.H.OR.ALPHA.GE.NP)  THEN 

F=1.0 

GOTO  29 
ELSE 

F=FUNC5(ALPHA,DELTA(H-K),EP) 
END  IF 

29 

XF=FUNC4(DSUM,U, ALPHA) 
IF  (U.LE.F)  THEN 

DT(I)=XF 

T(I)=T(I-1)+DT(I) 

c 

c 

IF  (T(I)-TIME(NTMAX2))  25,77,77 

IF  A  FAILURE  OCCURS  IN  THE  INTERVAL  WHERE  THERE  IS 

c 

ONLY  ONE  PLANT  OPERATING  THEN  THE  SYSTEM  CANNOT  FAIL 

c 

c 

25 

UNTIL  AT  LEAST  ONE  PLANT  RETURNS  TO  SERVICE. 

IF  (ALPHA. EQ.1)  THEN 

1=1+1 

DT(I)=BETA+T(I-NP)-T(I-1) 

T(I)=DT(I)+T(I-1) 

IF  (T(I)-TIME(NTMAX2))  26,77,77 

26 

IFLAG=1 
GOTO  5 
ELSE 
END  IF 
GOTO  5 
ELSE 
END  IF 

c 

75 

CONTINUE 

c 

THE 

LENGTH  OF  THE  INTERVAL  HAS  NOW  BEEN  DETERMINED  FOR  THE  CASE 
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C     DT(I-1)[BETA.  CONTINUE  IN  THE  PRESENT  SIMULATION. 

C     

GOTO  5 
77  RETURN 
END 
C 
C 
C 
C 

c  tittmtttntttttt* 

C  0SUBROUTINE  RESULT* 

C  M#f#f#####M####M 

SUBROUTINE  RESULT 

COMMON  T(100),SUM2(100),RATI02(100),TIME(100),TUSED(100),DOWN1(100 
1 ) ,DOWN2( 100) ,S1 (100) .RATI0 1 ( 100) ,SUM1 ( 100) .PFTR.BETA.I.NP.MAXIT.MO 
1DEL.NTMAX.NTMAX2 

C0MM0N/ERR0R/R1 ( 100) ,R2( 100) ,R3( 100) ,RU( 100) ,SA1 ( 100) ,SA2( 100) ,SA3 
1(100), SA4(100) 
COMMON/REPAIR/DEV 

C     

C     CALCULATE  THE  PROBABILITY  OF  PLAN  B  AT  THE  FAILURE  TIMES  FOR 
C     METHOD  1  OR  THE  FRACTION  DOWN  TIME  FOR  METHOD  2. 

C     

DO  20  J=1,NTMAX2,1 

TIME(J)=TIME(J)/BETA 
IF  (MODEL. EQ. 3)  GOTO  3 

C         

C         THE  RESULTS  FOR  MODEL  1  AND  2  ARE  CALCULATED  SIMULTANEOUSLY. 

C 

R ATIO 1 ( J ) = (MAXIT-SUM 1 ( J ) ) /MAXIT 

RSIGMA=SQRT(RATI01(J)»(1-RATI01(J))/MAXIT)«DEV 

R1(J)=RATI01(J)+RSIGMA 

R2 ( J ) =RATI0 1 ( J ) -RSIGMA 

SS= (MAXIT*BETA+S 1 ( J ) ) /MAXIT/BETA 

SSIGMA=SQRT ( SS» ( 1 -SS ) /MAXIT ) *DEV 

DOWN1(J)=1.0O-SS 

SA 1 ( J ) = 1 . OO-SS+SSIGMA 

SA2(J)=1.00-SS-SSIGMA 

RATI02(J)=(MAXIT-SUM2(J))/MAXIT 

RSIGMA=SQRT(RATI02(J)*(1-RATI02(J))/MAXm«DEV 

R3(J)=RATI02(J)+RSIGMA 

R1(J)=RATI02(J)-RSIGMA 

SS=RATI02(J) 

D0WN2(J)=1-SS 

SSIGMA=RSIGMA 

SA3(J)=1.00-SS+SSIGMA 

SA4(J)=1.00-SS-SSIGMA 

GOTO  20 

c         

C         THE  RESULTS  FOR  MODEL  3  ARE  CALCULATED. 

C         

3     RATI01(J)=(MAXIT-SUM1(J))/MAXIT 

RSIGMA=SQRT(RATI01(J)»(1-RATI01(J))/MAXIT)*DEV 
R1(J)=RATI01(J)+RSIGMA 
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R2(J)=RATI01(J)-RSIGMA 

SS= (MAXIT»BETA+S 1 ( J ) ) /MAXIT/BETA 

SSIGMA=SQRT(SS*(1-SS)/MAXIT)*DEV 

DOWN1(J)=1.00-SS 

SA 1 ( J ) = 1 . OO-SS+SSIGMA 

SA2(J)=1.00-SS-SSIGMA 
20  CONTINUE 
RETURN 
END 


****************** 

#SUBROUTINE  PRINT* 
****************** 

SUBROUTINE  PRINT 

COMMON  T(100),SUM2(100),RATI02(100),TIME(100),TUSED(100),DOWN1(100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

C0MM0N/ERR0R/R1 ( 100) ,R2( 100) ,R3( 100) ,R4( 100) ,SA1 ( 100) ,SA2( 100) ,SA3 
1(100),  SA4O00) 

COMMON/MAXT/TMAX ,MAXTAU ,NREAD ,NIBTR 
COMMON/REPAIR/DEV 
NFLAG=0 
10  WRITE  (6,20) 
20  FORMAT  (  '  1 ' , '  ##########OUTPUT  FROM  SIMULATION  FORTS AK##########> ) 

WRITE  (6,30)  BETA 
30  FORMAT  CO', 'REPAIR/REPLACEMENT  TIME  =  *,F10.6) 

WRITE  (6,50)  NP 
50  FORMAT  ('  ',' NUMBER  OF  PLANTS  IN  THE  POOL  =',I5) 

WRITE  (6,60)  MAXTAU 
60  FORMAT  ('  ', 'MAXIMUM  TIME  VALUE  (IN  UNITS  OF  BETA)  =',I5) 

WRITE  (6,70)  NIBTR 
70  FORMAT  (•  ', 'NUMBER  OF  INTERVALS  BEFORE  BETA  =',I4) 

WRITE  (6,80)  NTMAX2 
80  FORMAT  (•  ', 'TOTAL  NUMBER  OF  INTERVALS  =',I4) 

WRITE  (6,90)  MAXIT 
90  FORMAT  ('  ',' TOTAL  NUMBER  SIMULATIONS  RUN  =',I10) 
WRITE  (6,100)  MODEL 
100  FORMAT  (•  ', 'MODEL', 12,'  WAS  USED') 

WRITE  (6,110)  DEV,DEV,DEV,DEV 
110  FORMAT  CO','  TIME  PROB. AVAIL  +',F4.2,'  SIGMA  -',F4.2,'  SIGMA  • 
1,'FRAC.  DOWN  +',F4.2,'  SIGMA  -',F4.2,'  SIGMA') 
WRITE  (6,120) 

120  FORMAT  ('  ','  • 

1 , ' . ) 

DO  140  J=1,NTMAX2,1 
IF  (MOD EL. EQ. 3)  GOTO  3 
IF  (MODEL. EQ. 2)  GOTO  2 

WRITE  (6,130)  TIME(J),RATI01(J),R1(J),R2(J),D0WN1(J),SA1(J),SA2(J) 
130  FORMAT (■  ',F6.3,'    ' ,F7. V    ',F9.6,'    ',F9.6,'    ',F8.6,'   ' 
1.'   ',F9.6,'    ',F9.6) 
GOTO  140 
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2  WRITE  (6,130)  TIME(J),RATI02(J),R3(J),R4(J),DOWN2(J),SA3(J),SA4(J) 
GOTO  140 

3  WRITE  (6,130)  TIME(J),RATI01(J),R1(J),R2(J),DOWN1(J),SA1(J),SA2(J) 
140  CONTINUE 

IF  (NFLAG.EQ.1)  GOTO  150 
IF  (MODEL.EQ.1)  THEN 

M0DEL=2 

NFLAG= 1 

GOTO  10 
ELSE 
END  IF 
IF  (MODEL. EQ. 2)  THEN 

M0DEL=1 

NFLAG=1 

GOTO  10 
ELSE 
END  IF 
150  RETURN 
END 


******************* 

^SUBROUTINE  RANDOM* 
******************* 

SUBROUTINE  ROUTINE  THAT  OUTPUTS  A  RANDOM  NUMBER  FROM  0  TO  1  FOR  A 
GIVEN  INPUT  SEED  VALUE. 

SUBROUTINE  RANDOM 

COMMON/RNG/DSEED,U 

DOUBLE  PRECISION  DSEED.U 

DSEED=DMOD( 1 6807. D0*DSEED, 2 147483647. DO) 

U=DSEED/21 47483648. DO 

RETURN 

END 


******************* 

#SUBROUTINE  CHECK 1# 
******************* 

THE  SUBROUTINE  CHECK  DETERMINES  IF  THE  PRESENT  VALUE  OF  TIME 
EXCEEDS  ANY  OF  THE  TIME  VALUES  AT  WHICH  THE  K-TH  PLANT  FAIL- 
URE OCCURS. 

SUBROUTINE  CHECK  1 

COMMON  T(100),SUM2(100),RATI02(100),TIME(100),TUSED(100),DOWN1(100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 

M=1 
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DO  90  J=1,I,1 

DO  80  L=M,NTMAX2,1 

IF  (T(J)-TIME(D)  85,10,10 
10         IF  (TUSED(J))  70,70,30 
30         IF  (TIME(D-TUSED(J)-BETA)  20,70,70 
20         SUM2(L)=SUM2(L)+1.0 
70         IF  (J.EQ.1)  GOTO  80 
40         TNEW=TIME(L)-T(J-1)-BETA 

IF  (TNEW)  60,80,80 
60         SUM1(L)=SUM1(L)+1.0 

S1(L)=S1(L)+TNEW 
80     CONTINUE 
85     M=L 
90  CONTINUE 
RETURN 
END 
C 
C 
C 
C 

C  ******************* 

C  #SUBR0UTINE  CHECK3* 

C  ******************* 

C     

C  THE  SUBROUTINE  CHECK  DETERMINES  IF  THE  PRESENT  VALUE  OF  TIME 
C  EXCEEDS  ANY  OF  THE  TIME  VALUES  AT  WHICH  THE  K-TH  PLANT  FAIL- 
C     URE  OCCURS. 

C     

SUBROUTINE  CHECK 3 

COMMON  T( 100) ,SUM2( 100) ,RATI02( 100) ,TIME( 100) ,TUSED( 100) ,DOWN1 ( 100 
1),DOWN2(100),S1(100),RATI01(100),SUM1(100),PFTR,BETA,I,NP,MAXIT,MO 
1DEL,NTMAX,NTMAX2 
M=1 
DO  90  J=1,I,1 

DO   80  L=M,NTMAX2,1 

IF    (T(J)-TIME(D)    85,10,10 
10  IF    (J.EQ.1)   GOTO   80 

TNEW=TIME(L)-T(J-1)-BETA 
IF   (TNEW)    60,80,80 
60  SUM1(L)=SUM1(L)+1.0 

S1(L)=S1(L)+TNEW 
80  CONTINUE 

85  M=L 

90  CONTINUE 
RETURN 
END 
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APPENDIX  B:   KSUSPARE  Computer  Code 

The  computer  code  KSUSPARE  was  developed  to  determine  the  present 
worth  of  Plan  A  and  Plan  B  costs.   The  computer  code,  written  in 
FORTRAN,  consists  of  a  main  program  and  13  subroutines. 

Input  Data 

The  main  program  calls  subroutines  to  read  in  and  print  out  the 
problem  data  and  to  perform  the  economic  analysis.   The  program  is 
currently  set  up  to  estimate  the  Plan  A  and  Plan  B  costs,  the  expected 
savings,  and  the  benefit-cost  ratio  for  one  plant  in  the  spare  pool  over 
the  spare  pool  lifetime  at  variable  time  steps. 

LAMBDA  =  component  failure  rate. 

TREPL  =  repair/replacement  time. 

SPRMW  =  spare  component  capacity  (MW) . 

N  =  number  of  plants  in  the  pool. 

NI  =  number  of  intervals  used  for  numerical  integration  in  all  but 
the  subroutine  GOT  (i.e.,  the  subroutine  to  estimate  g(t)). 

NIG  =  number  of  intervals  used  for  numerical  integration  in  the 
subroutine  GOT. 

METHOD  =  the  values  of  g(t)  and  f   (t)  are  constant  when  a  value  of 
1  is  read  in,  otherwise  these  values  are  determined  using 
numerical  integration. 

IT  =  the  maximum  evaluation  time  as  an  integer  value. 

TI  =  the  maximum  evaluation  time  as  a  real  value. 

TLOOP  =  the  spare  pool  lifetime. 

QSPR  =  present  worth  cost  of  the  spare. 

QSSV  ■  present  worth  salvage  value  of  the  spare. 
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QSTOR  ■  present  worth  storage  cost  per  year  for  the  spare. 
MODEL  =  selects  the  method  by  which  the  fraction  downtime  for  each 
time  interval  is  assigned.   An  input  value  of  4  causes  the 
fraction  downtime  to  be  calculated  numerically  using 
renewal  theory;  otherwise,  the  fraction  downtime  is  read 
in  by  data  statements. 
The  remaining  input  parameters  are  plant-specific  costs  and  data. 

Subroutines 


The  subroutine  PWORTH  calls  other  subroutines  to  calculate  the  Plan 
A  and  Plan  B  present  worth  costs  for  each  plant.   The  subroutine  C0ST1 
calculates  the  failure-dependent  variable  and  fixed  costs  and  the 
variable  costs  during  plant  operation.  The  subroutine  C0ST2  calculates 
the  annual  fixed  costs.   The  subroutine  GOT  calculates  the  numerical 
solution  to  the  renewal  equation  and  returns  the  value  of  g(t)  for  each 
time  interval.   The  subroutine  DCOUNT  and  ESCLAT  calculate  the  discount 
and  escalation  rates  for  each  time  interval.   The  subroutine  CPC1TY 
initializes  the  interval  values  for  plant  capacity  factors  and  fraction 
downtime  for  each  time  interval.   The  fraction  downtime  is  calculated 
using  renewal  theory  by  subroutine  DTIME  if  Model  4  is  selected. 

The  escalation  and  discount  rates  are  initialized  in  their 
respective  subroutines  by  DATA  statements.   Thus,  for  each  simulation 
model  the  data  statements  must  be  modified  and  a  separate  version 
compiled.   The  remaining  data  are  read  in  from  an  external  data  file. 

Code  Listing 

The  following  is  a  listing  of  the  KSUSPARE  computer  code. 
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$JOB         VAS,TIME=(2,40),PAGES=25 

n***********»x*****ll*********KS[jSPAi;  £•»******  ******  »****»****»«**»****** 

c«  * 

C»  THIS  PROGRAM  ESTIMATES  THE  TOTAL  PRESENT  WORTH  COSTS  FOR  PLAN  A  * 
C»  AND  PLAN  B.  THE  PROGRAM  USES  THE  TRAPEZOID  RULE  WHEN  NUMERICAL  » 
C*    INTEGRATION  IS  PERFORMED.  » 

C»  » 

r  **»***********»**»»****»**»»**♦  *»**♦«*♦**»*«**«****»»****»*******»»»*»* 

C  *MAIN  PROGRAM* 

COMMON  DTF(1000),CAP(1000),ENA,ENB,OMA,OMB,REV,SRC,FOR,TREV,TLOOP, 
1TOTALA,TOTALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QSTOR,QFSV,QPCSV,QSPR, DELTA, SPRMW.MW.TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 


REAL  LAMBDA, MW 

CALL  READ1 

CALL  PRINT 1 

CALL  PWORTH 

STOP 

END 

c 

c 

c 

c 

c 

c 

c 

•SUBROUTINE  READ1* 

SUBROUTINE  READ1 

COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV.TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,C0MP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
IMP 

REAL  LAMBDA, MW 

READ  (5,10)  LAMBDA, TREPL,TTEMP,SPRMW 
10  FORMAT  (/,F10.3,2F10.6,F10.3) 

'READ  (5,20)  N,NI,NIG,METHOD,MODEL 
20  FORMAT  (/,5I5) 

READ  (5,30)  QSPR,QSSV,QSTOR 
30  FORMAT  (/.3F10.3) 

READ  (5,40)  IT,TI,STEP,TLOOP 
10  FORMAT  C/,I5,3F10.3) 
RETURN 
END 
C 
C 
C 
C 
C  *•••****«***»**•••• 

C  'SUBROUTINE  PRINT  1* 
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SUBROUTINE  PR INT 1 

COMMON  DTF ( 1 000 ), CAP (1 000 ),ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV, TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR,LAMBDA,TREPL,QEN,Q0M,QC0MP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MOD EL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA, MW 
WRITE  (6,5) 
5  FORMAT  ( ' 1 ' , •###################KSUSPARE########################' ) 
WRITE  (6,10)  LAMBDA 
10  FORMAT  ('  ', 'COMPONENT  FAILURE  RATE  (YR-1)  =',F7.3) 

WRITE  (6,20)  TREPL 
20  FORMAT  (■  ' , 'NEW  COMPONENT  REPLACEMENT  TIME  (YR)  =',F10.6) 

WRITE  (6,30)  IT 
30  FORMAT  (•  '.'MAXIMUM  EVALUATION  TIME  (YR)  =',I5) 

WRITE  (6, UO)  STEP 
40  FORMAT  ('  ', 'EVALUATION  TIME  STEP  (YR)  =',F7.3) 

WRITE  (6,50)  TLOOP 
50  FORMAT  ('  ', 'SPARE  POOL  LIFETIME  (YR)  =',F10.3) 

WRITE  (6,60)  N 
60  FORMAT  ('  ', 'NUMBER  OF  PLANTS  IN  THE  SPARE  POOL  =',I5) 

WRITE  (6,70)  NI 
70  FORMAT  (•  ', 'NUMBER  OF  INTEGRATION  INTERVALS  (EXCEPT  GOT)  =',I5) 

WRITE  (6,80)  NIG 
80  FORMAT  ('  ', 'NUMBER  OF  INTEGRATION  INTERVALS  FOR  SUB.  GOT  =',I5) 

WRITE  (6,90)  METHOD, MODEL 
90  FORMAT  (■  ','G(T)  METHOD  NUMBER  ',11,'  AND  MODEL', 12,'  WERE  USED') 
WRITE  (6,100)  SPRMW 
100  FORMAT  (•  ', 'SPARE  COMPONENT  CAPACITY  (MW)  =',F10.3) 

WRITE  (6,110)  QSPR 
110  FORMAT  (•  ', 'SPARE  COMPONENT  COST  (K$)  =',F10.3) 

WRITE  (6,120)  QSSV 
120  FORMAT  ('  ', 'SPARE  COMPONENT  SALVAGE  VALUE  (K$)  =',F10.3) 
RETURN 
END 
C 
C 
C 
C 

Q  *»******»»***«»»*»* 

C  "SUBROUTINE  PWORTH* 

C  »»»»»»»»»»»»»»*»«»» 

SUBROUTINE  PWORTH 

COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV, TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL, QEN,QOM,QCOMP, QSSV, QTEMP, 
1QRM,QST0R,QFSV,QPCSV, QSPR, DELTA, SPRMW, MW.TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 
REAL   LAMBDA, MW 

C  

C  THE  INTERVAL  WIDTH  FOR  THE  NUMERICAL   INTEGRATION  IS  CALCULATED. 

C  

Y=FLOAT(NI) 
DELTA=1/Y 
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C     

C     FBACTION  OF  YEARLY  FIXED  COSTS  PAID  BY  EACH  PLANT  IS  CALCULATED. 
C     

YN=FLOAT(N) 

FRAC=1/YN 

ICOUNT=0 

NP=1 

C     

C     SPARE  COMPONENT  SALVAGE  VALUE  IS  ESTIMATED. 
C     


CALL  READ2 

CALL  PRINT2 

CALL  CPCITY 

CALL  REVNUE 

CALL  ESCLAT  (TLOOP.EF1 ,EF2,EF3,EF4,EF5,EF6,EF7,EF8,EF9) 

C     

C     ALL  DATA  HAS  BEEN  INITIALIZED. 
C     


ICOUNT=1 

CALL  DCOUNT  (TLOOP.DF) 

SSV=DF  »FR  AC  *Q  SS V*R  E V*EF9 

C 

C     SPARE  COMPONENT  PURCHASE  COST  IS  ESTIMATED . 
C     

SPR=QSPR»FRAC*REV 


C  THE  PLANT  DEPENDENT  SUBROUTINES  ARE  NOW  CALLED   AFTER  THE  PLANT  TO 

C  BE  ANAYLZED  IS  ESTABLISHED. 


T=0.0 
C  

C     THE  OUTSIDE  LOOP  IS  SETUP  HERE. 
C     

T=STEP 
10  CALL  C0ST1 

CALL  C0ST2 

CALL  RESULT 

CALL  PRINT3 

T=T+STEP 

IF  (T.GT.TI)  GOTO  20 

GOTO  10 
20  CONTINUE 

RETURN 

END 
C 
C 
C 
C 
C  **•>****•***•**•*• 

C  "SUBROUTINE  READ2* 

C  *«*•***««*«•••**•• 

SUBROUTINE  READ2 

COMMON  DTF( 1000), CAP( 1000), ENA.ENB.OMA, OMB, REV, SRC, FOR, TREV.TLOOP, 
1TOTALA,T0TALB,SAV,BCR,T,FCR,LAMBDA,TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
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1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA, MW 

READ    (5,10)   QEN,QOM,QCOMP,QPCSV,QFSV,QTEMP,QRM 
10  FORMAT    (/,7F10.3) 

READ    (5,20)   MW,FCR,FOR,SRC,TREV 
20  FORMAT   (/,5F10.3) 

RETURN 

END 


•SUBROUTINE  PRINT2* 

SUBROUTINE  PRINT2 

COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA,OMB, REV, SRC, FOR, TREV.TLOOP, 
1TOTALA,TOTALB,SAV,BCR,T,FCR,LAMBDA,TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QSTOR,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA.MW 
WRITE   (6,10)    NP 

10  FORMAT    ('-',' INPUT/OUTPUT  FOR   PLANT  NUMBER    ',13,' ') 

WRITE   (6,20) 

20  FORMAT    ('    ', 'ENERGY O&M COMP.-COMP.    SAL. -FAILED  COMP.    SAL' 

1 , '  —REPLACEMENT-RES. MARGIN' ) 

WRITE    (6,30)   QEN,QOM,QCOMP,QPCSV,QFSV,QTEMP,QRM 
30  FORMAT    (•    ',9F10.3) 

WRITE    (6,40)   MW.FCR.TREV 
40  FORMAT    ('    ', 'CAPACITY    (MW)    =',F10.3,'      FCR   =',F7.3,'      YEARS    (REV)' 
1,'    =',F7.3) 
WRITE   (6,50)    FOR, SRC 
50  FORMAT    (*    ', 'FORCED  OUTAGE  RATE   =',F7.3,'      SYSTEM  RISK  CHARACTERS 
1,'STIC   =',F10.3) 
WRITE   (6,60) 
60  FORMAT    ('-','    TIME     TOT ALA     TOTALB  SAV        BCR  ENERGY  ■ 

1,'  O&M  RESERVE  COMP        FSV  USV        SSV     TEMP      ' 

1,' SPARE        STOR') 
WRITE   (6,70) 

70  FORMAT    (*    ',' ' 

1 ,  ■ . 

1,. •) 

RETURN 
END 


I***************** 
•SUBROUTINE  COST1* 
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SUBROUTINE  COST1 

COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV,TLOOP, 
1 TOTALA , TOT ALB , SAV , BCR , T , FCR , LAMBDA , TREPL , QEN , QOM , QCOMP , QSSV , QTEMP , 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,C0MP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT .METHOD, MODEL, NP,N,IC0UNT,TE 
1MP 

DIMENSION  VCA( 1000) ,VCB( 1000) ,FCA( 1000) ,FCB( 1000) ,FCC( 1000) ,VCC( 10 
100), VCD( 1000) 
REAL  LAMBDA, MW 

C     

C     THE  CONSTANT  MULTIPLIERS  ARE  CALCULATED  FIRST. 
C     


C1A=TREPL»QEN*LAMBDA 

C1B=(MW-SPRMW) 

IF   (C1B.LT. 0.0)    C1B=0.0 

C2=TREPL*Q0M 

IF  (METH0D.EQ.1)  THEN 

C3=1/( 1+LAMBDA»TREPL) 
ELSE 

C3=1.0 
END  IF 

C     

C     THE  VALUES  FOR  THE  FUNCTION  TO  BE  NUMERICALY  INTEGRATED  ARE  FIRST  EVAL- 
C     UATED  AT  THE  TIME  VALUES  CALCULATED  ABOVE. 

C     

ID=NI*T 
ID2=ID+1 
VCA(1)=CAP(1) 
VCB(1)=0.0 
DO    10  K=2,ID2,1 
J=K-1 
X=J«DELTA 
CALL  GOT(X.GTX) 
CALL  ESCLAT    (X,EF1 ,EF2,EF3,EF4,EF5,EF6,EF7,EF8,EF9) 

C  

C         THE  LAST  PIECE  OF  DATA  HAS  BEEN  INITIALIZED. 
C         


IC0UNT=1 

CALL  DCOUNT  (X,DF) 

VCA(K)=GTX«DF*EF1*CAP(J) 

VCB(K)=VCA(K)*DTF(J) 

VCC(K)=DF«EF2*GTX 

VCD(K)=DF»EF2*GTX*(1-DTF(J)) 

FCA(K)=GTX»EF3*DF 

FCB(K)=GTX«EF4«DF 

FCC(K)=GTX*EF5*DF 
10  CONTINUE 

C     

C     THE  NUMERICAL  INTEGRATION  IS  NOW  PERFORMED. 
C     


SUM1=0.0 
SUM2=0.0 
SUM3=0.0 
SUM5=0.0 
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SUM6=0.0 
SUM7=0.0 
SUM8=0.0 
SUM9=0.0 
DO   20  K=2,ID,1 

SUM1=SUM1+VCA(K) 
SUM2=SUM2+VCB(K) 
SUM5=SUM5+FCA(K) 
SUM6=SUM6+FCB(K) 
SUM7=SUM7+FCC(K) 
SUM8=SUM8+VCC(K) 
SUM9=SUM9+VCD(K) 
20  CONTINUE 

C  

C     THE  INITIAL  VALUES  FOR  THE  F(0)  TERMS  USED  IN  THE  NUMERICAL  INTE- 
C     GRATION  ARE  ESTABLISHED.  IF  THE  CONSTANT  METHOD  IS  USED  THE  F(0) 
C     TERM  IS  SET  EQUAL  TO  SOME  CONSTANT  VALUE. 

C     

IF  (METH0D.EQ.1)  THEN 

GT=1/(1+LAMBDA*TREPL) 
VCA(1)=CAP(1) 
VCB(1)=DTF(1) 
VCC(1)=GT 
VCD(1)=GT 
FCA(1)=GT 
FCB(1)=GT 
FCC(1)=GT 
ELSE 

VCA(1)=CAP(1) 
VCB(1)=0.0 
VCC(1)=1.0 
VCD(1)=1.0 
FCA(1)=1.0 
FCB(1)=1.0 
FCC(1)=1.0 
END   IF 

C  

C            ENERGY   REPLACEMENT  COST  IS  ESTIMATED. 
C  


AA=DELTA/2*(VCA(1)+VCA(ID2)+2«SUM1) 

BB=DELTA/2«(VCB(1)+VCB(ID2)+2*SUM2) 

ENA=MW»C1A*AA 

ENB=C 1 A« ( C 1 B«AA+SPRMW*BB) 

C  

C     OPERATION  AND  MAINTENANCE  COSTS  ARE  ESTIMATED. 
C     

CC=DELTA/2»( VCC( 1 )+VCC(ID2)+2»SUM8) 

DD=DELTA/2*(VCD(1)+VCD(ID2)+2»SUM9) 

OMA=QOM«CC 

OMB=OMA+LAMBDA*TREPL*DD*QOM 

C  

C            NEW  COMPONENT   PURCHASE  COST   IS  ESTIMATED. 
C  

CC=DELTA/2*(FCA( 1 )+FCA(ID2)+2*SUM5) 
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COMP=QCOMP«CC«REV«LAMBDA 
C  

C            SALVAGE  OF  THE  FAILED  COMPONENT   IS  ESTIMATED. 
C  


DD=DELTA/2«(FCB(1)+FCB(ID2)+2»SUM6) 
FSV=QFSV«DD*REV*LAMBDA 

C  

C  TEMPORARY   REPLACEMENT  COST  IS  ESTIMATED. 

C  

IF    (MW-SPRMW.LE.0.0)    THEN 

TEMP=0.0 
ELSE 

EE=DELTA/2»(FCC( 1 )+FCC(ID2)+2»SUM7) 
TEMP=LAMBDA»EE*QTEMP 
END   IF 

C  

C     PLANT  COMPONENT  SALVAGE  VALUE  IS  ESTIMATED. 
C     


CALL  DCOUNT  (T,DF) 

CALL  GOT  (T.GTX) 

CALL  ESCLAT  (T,EF1 ,EF2,EF3,EF4,EF5,EF6,EF7,EF8,EF9) 

PCSV=QPCSV*EF8«DF*REV*GTX 

RETURN 

END 
C 
C 
C 
C 
C  •»*»•»»»*»*»**»««» 

C  'SUBROUTINE  COST2* 

C  ****«***«*«***»*»» 

SUBROUTINE  C0ST2 

COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA,OMB, REV, SRC, FOR, TREV,TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC.STOR ,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG,  IT  ,-METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

DIMENSION  XX( 1000), YY( 1000) 
REAL  LAMBDA.MW 

C 

C     RESERVE  CAPACITY  REQUIRED  CALCULATION  IS  FIRST.  THE  CONSTANTS  USED 
C     IN  DETERMINING  THE  DISCOUNT/ESCALATION  FACTOR  AND  IN  NUMERICALY 
C     EVALUATING  THE  AVERAGE  DOWN  TIME  ARE  ESTIMATED. 

C     

IF  (T.LE.1.0)  THEN 

CALL  DCOUNT  (1.0,DF) 

CALL  ESCLAT  ( 1 .0,EF1 ,EF2,EF3,EF4,EF5,EF6,EF7,EF8,EF9) 
C1=DF*EF6 
C2=DF«EF7 

IF  CMETH0D.EQ.1)  THEN 
FBAR=DTF(1) 
GOTO  UO 
ELSE 
END  IF 
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BB=T*NI 
NT1=INT(BB)+1 
NB1=2 
NB2=NB1 
XX(1)=1.0 
YY(1)=1.0 
L=1 
ELSE 

X=T-INT(T) 

IF   (X.EQ.O.O)    THEN 

M=INT(T) 
ELSE 

M=INT(T)+1 
END   IF 

NB1=(M-1)*NI+1 
NT1=M»NI+1 
NB2=NB1+1 
L=NB1 
SUMUO.O 
SUM2=0.0 

C  

C         DISCOUNT/ESCALATION  FACTOR  FOR  RESERVE  MARGIN. 

C         

DO    10  K=1,M,1 
X=FLOAT(K) 
CALL  DCOUNT   (X,DF) 

CALL  ESCLAT    (X.EF1  ,EF2,EF3,EF1»,EF5,EF6,EF7,EF8,EF9) 
SUM1=SUM1+DF*EF6 
SUM2=SUM2+DF*EF7 
10  CONTINUE 

C1=SUM1 
C2=SUM2 
END   IF 

C  

C     THE  AVERAGE  DOWN  TIME  FRACTION  FOR  YEARLY  INTERVALS  IS  NUMERICALY 
C     EVALUATED. 

C     

DO  20  K=NB1,NT1,1 
J=K-1 
S=K*DELTA 
CALL  GOT(S,GTS) 
XX(K)=GTS*DTF(J) 
YY(K)=GTS 
20  CONTINUE 
NT2=NT1-1 
SUM5=0.0 
SUM6=0.0 

DO   30  K=NB2,NT2,1 
SUM5=SUM5+XX(K) 
SUM6=SUM6+YY(K) 
30  CONTINUE 

FBAR=(XX(L)+XX(NT1)+2«SUM5)/CYY(L)+YY(NT1)+2*SUM6) 

C  

C     THE  RESERVE  MARGIN  DUE  TO  THE  COMPONENT  AND  ASSOCIATED  COST  IS 
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C  ESTIMATED. 

C  

40  RCA=LAMBDA*TREPL«(1-FOR) 

AA=EXP(MW/SRC)-1 

RATIO=(MW-SPRMW)/MW 

IF    (RATIO. LT. O.O)    RATIO=0.0 

RMA=SRC»AL0G(1+(RCA«AA/(1+F0R*AA)))»FCR«QRM*C1 

RCBzLAMBD A»TREPL« ( FBAR+( 1 -FBAR ) "RATIO ) » ( 1 -FOR ) 

RMB=SRC*ALOG ( 1 +( RCB*AA/ ( 1 +FOR*AA ) ) ) »FCR«QRM»C 1 

C  

C     THE  SPARE  COMPONENT  STORAGE  COST  IS  ESTIMATED. 
C     

ST0R=QST0R*FRAC*C2 

RETURN 

END 
C 
C 
C 
C 
C  »•»»»»»*»*«»»»•»»»» 

C  'SUBROUTINE  RESULT* 

C  »*»*»«*****♦»**»*** 

SUBROUTINE  RESULT 

COMMON  DTF( 1000), CAP( 1000), ENA.ENB.OMA, OMB, REV, SRC, FOR, TREV,TLOOP, 
1T0T ALA, TOTALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT .METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA.MH 
C 
C 
C 

C 
C 
C 

T0TALB=ENB+OMB+C0MP-FSV+TEMP+RMB+ST0R-PCSV-SSV+SPR 

C  

C     THE  SAVINGS  IS  THE  DIFFERENCE  BETWEEN  THE  TWO  PLANS. 


DETERMINE 

THE 

ESTIMATE 

OF 

THE 

TOTAL 

PRESENT 

WORTH 

COST 

FOR 

PLAN 

A. 

TOTALA=ENA+OMA+COMP-FSV+RMA-PCSV 

DETERMINE 

THE 

ESTIMATE 

OF 

THE 

TOTAL 

PRESENT 

WORTH 

COST 

FOR 

PLAN 

B. 

c     

SAV=TOTALA-TOTALB 
C     

C     THE  FINAL  STEP  IS  TO  CALCULATE  THE  BENEFIT-COST  RATIO. 
C     

BCR=1+(SAV/SPR) 

RETURN 

END 
C 
C 

c 
c 

C  »»»*«»»»*»»»»«»»»«» 

C  "SUBROUTINE  PRINT3* 
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SUBROUTINE   PRINT3 

COMMON  DTF( 1000), CAP( 1000), ENA, ENB, OMA, OMB, REV, SRC, FOR, TREV.TLOOP, 
1T0TALA, TOT ALB, SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA, MW 

WRITE   ( 6 , 30 )    T , TOTALA , TOTALB , SAV , BCR , ENA , ENB , OMA , OMB , RMA , RMB , COMP , 
1FSV,PCSV,SSV,TEMP,SPR,ST0R 
30  FORMAT   ('    ',F5.2,3F8.2,f6.2,F8.2,F7.2,3F8.2,2F7.2,f6.2,F7.2,2F6.2, 
1F7.2,F6.2,F8.2,F6.2) 
RETURN 
END 
C 
C 
C 
C 

Q  ***«****«*«****» 

C  'SUBROUTINE  GOT* 

SUBROUTINE  G0T(Y,GTY) 

COMMON  DTF( 1000), CAP( 1000), ENA, ENB, OMA, OMB, REV, SRC, FOR, TREV,TLOOP, 
1T0TALA, TOTALB, SAV, BCR, T.FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW.MH.TI, STEP, FRAC,ST0R, COMP, SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI,NIG,IT,METH0D,M0DEL,NP,N,IC0UNT,TE 
1MP 

DIMENSION  GTC1000) 
REAL  LAMBDA, MW 

C  

C     THE  VALUE  OF  G(T)  CAN  EITHER  BE  CONSTANT  OR  ESTIMATED  BY  NUMERICAL 
C     INTEGRATION.  IF  THE  METHOD  IS  1  THEN  G(T)  IS  CONSTANT. 

C     

IF  (METH0D.EQ.1)  THEN 

GTY=1/( 1+LAMBDA*TREPL) 
GOTO  50 
ELSE 
END  IF 
IF  (Y.LE.TREPL)  THEN 

GTY=EXP(-LAMBDA*Y) 
ELSE 

C         

C         CONSTANTS  USED  IN  THE  SUBROUTINE  ARE  CALCULATED  FIRST. 

C         

C1=EXP(TREPL) 

ALPHA=LAMBDA»EXP(LAMBDA«TREPL) 

DELTA2=TREPL/NI 

BETA=ALPHA»DELTA2/2 

ID3=Y/DELTA2+1 

NIUNI  +  1 

DO    10  1=1 ,NI1 , 1 

X=(I-1)«DELTA2 

C  

C            THE  FIRST  TREPL  INTERVAL  IS  ESTIMATED. 
C  
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GT(I)=EXP(-LAMBDA*X) 
10  CONTINUE 

NI2=NI1+1 
DO   HO  I=NI2,ID3,1 
X=(I-1)*DELTA2 

C  

C  THE  REMAINING  TIMES  ARE  ESTIMATED. 

C  

AA=EXP(-LAMBDA»X) 

SUM=0.0 

NEND=I-NI2 

IF  (NEND.EQ.O)  GOTO  30 

DO  20  J=1,NEND,1 

SUM=SUM+EXP(LAMBDA«(J)*DELTA2)*GT(J+1) 
20        CONTINUE 
30         BB=EXP(LAMBDA*(X-TREPL))*GT(I-NI) 

GT(I)=AA«(1+BETA*(GT(1)+BB+2*SUM)) 
40     CONTINUE 

GTY=GT(ID3) 
END  IF 
50  CONTINUE 
RETURN 
END 
C 
C 
C 
C 

C  "SUBROUTINE  DCOUNT* 

Q  •*****••***•***•••* 

SUBROUTINE  DCOUNT  (Y,DF) 

COMMON  DTFC 1000), CAP( 1000), ENA,ENB,OMA,0MB, REV, SRC, FOR, TREV,TLO0P, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI,NIG,IT,METH0D,M0DEL,NP,N,IC0UNT,TE 
1MP 

DIMENSION  DR(100) 
REAL  LAMBDA, MW 

C  

C     DISCOUNT  RATES  MUST  BE  ASSIGNED  FOR  THE  DIFFERENT  PLANTS.  IF  THE 
C     VALUES  HAVE  ALREADY  BEEN  INITIALIZED  THEN  THE  DATA  STATEMENT  IS 
C     BYPASSED. 

C     

IF  (IC0UNT.EQ.1)  GOTO  10 
DATA  DR/ 100*0. 13/ 
10  M=INT(Y) 

IF  (Y.LE.1.0)  M=0 
DF=(1+DR(M+1))«*(-Y) 
RETURN 
END 
C 
C 
C 

c 
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Q  MM*************** 

C  'SUBROUTINE  ESCLAT* 

Q  ******************* 

SUBROUTINE  ESCLAT(X,XEF1 ,XEF2,XEF3,XEF4,XEF5,XEF6,XEF7,XEF8,XEF9) 
COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV.TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

DIMENSION  ER1 ( 100) ,ER2( 100) ,ER3( 100) ,ER4( 100) ,ER5( 100) ,ER6( 100) ,ER 
17(100),ER8(100),ER9(100) 
REAL  LAMBDA, MW 
IF    (IC0UNT.EQ.1)   GOTO    10 

C  

C            ENERGY  RATE 
C  


DATA  ER 1/100*0.04/ 
C  

C     OPERATION  AND  MAINTENANCE  RATE 

C     

DATA  ER2/1 00*0.03/ 

C     

C     NEW  COMPONENT  RATE 
C     


DATA  ER3/5*0.05,95*0.07/ 

C     

C     FAILED  COMPONENT  SALVAGE  RATE 
C     

DATA  ER<t/ 100*0. 055/ 

C  

C            TEMPORARY   REPLACEMENT  RATE 
C  

DATA  ER5/ 100*0. 035/ 

C  

C     RESERVE  MARGIN  RATE 
C     


DATA  ER6/ 100*0. 05/ 
C     

C     SPARE  STORAGE  RATE 
C     

DATA  ER7/ 100*0. 03/ 
C 
C 
C 

C 
C 
C 

DATA  ER9/ 100*0. 055/ 

C  

C  IF  THE  TIME  TO  EVALUATE  THE  ESCALATION  RATE  IS  LESS  THAN  ONE  YEAR 

C  THEN  THE  ESCALATION  RATE  IS  BASED  ON  A  FRACTION  OF  THE  FIRST  YEAR. 

C  


USED 

COMPONENT 

SALVAGE 

VALUE 

RATE 

DATA 

ER8/ 100*0. 

055/ 

SPARE  COMPONENT  SALVAGE 

VALUE 

RATE. 

10  M=INT(X) 
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»**««*«*«******»»»* 


IF    (X.LE.1.0)   M=0 

XEF1=(1+ER1(M+1))*«X 

XEF2=(1+ER2(M+1))»«X 

XEF3=(1+ER3(M+1))**X 

XEF4=(1+ER4(M+1))**X 

XEF5=(1+ER5(M+1))»«X 

XEF6=(1+ER6(M+1))«»X 

XEF7=(1+ER7(M+1))»«X 

XEF8=(1+ER8(M+1))»*X 

XEF9=(1+ER9(M+1))«*X 

RETURN 

END 
C 
C 
C 
C 
C 

C  'SUBROUTINE  REVNUE* 

C  ******♦***»******«» 

SUBROUTINE  REVNUE 

COMMON  DTF( 1000), CAP( 1 000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV.TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SFRMW,MW,TI, STEP, FR AC, STOR.COMP.SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI,NIG,IT,METH0D,M0DEL,NP,N,IC0UNT,TE 
1MP 

REAL  LAMBDA, MW 

SUM2=0.0 

X=0.0 
10  X=X+1.0 

CALL  DCOUNT    (X,DF) 

SUM2=SUM2+DF 

IF    (X.LT.TREV)   GOTO    10 

REV=FCR«SUM2 

RETURN 

END 

C 

c 
c 

Q  »»»»»«»»»»»*»»«»»»» 

C  'SUBROUTINE  CPCITY* 

SUBROUTINE  CPCITY 

COMMON  DTFC 1000), CAP( 1000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV.TLOOP, 
1 TOTALA , TOT ALB , SAV , BCR , T , FCR , LAMBDA , TREPL , QEN , QOM , QCOMP , QSSV , QTEMP , 
1QRM,QST0R,QFSV,QPCSV,QSPR, DELTA, SPRMW,MW,TI, STEP, FRAC,STOR ,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA, MW 
DIMENSION  FDT  C 1 00 ) , FC AP ( 1 00 ) 

C     

C     ASSIGN  CAPACITY  FACTORS  AND  FRACTION  DOWN  TIMES  FOR  THE  DIFFERENT 
C     PLANTS.  IF  THE  VALUES  HAVE  ALREADY  BEEN  ASSIGNED  THEN  THE  DATA 
C     STATEMENT  IS  BYPASSED. 
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IF  (IC0UNT.EQ.1)  GOTO  10 

DATA  FCAP.FDT/ 100*0. 390, 0.04 18, 0.0429, 0.0430, 0.0427, 0.044, 0.043,0. 
10433, 0.0444, 0.0435, 0.0417, 0.0414, 0.0419, 0.043, 0.0397, 0.043,85*0.04 
13/ 
10  ID=IT«NI 

DO  20  1  =  1, ID,  1 
X=I*DELTA 
M=INT(X) 
IF  (MOD EL. EQ. 4)  THEN 

CALL  DTIME  (X,DTX) 

DTF(I)=DTX 

IF  (METH0D.EQ.1)  THEN 

CAP(I)=FCAP(1) 
ELSE 

R=M-IT 

IF  (R.EQ.O.O)  M=M-1 
CAP(I)=FCAP(M+1) 
END  IF 
ELSE 


IF  THE  THE  EVALUATION 

TIME  IS 

LESS 

THAN 

ONE 

YEAR 

THEN 

THE 

CAPACITY 

FACTOR 

AND  FRACTION 

DOWN 

TIME 

ARE 

BASED 

ON  A 

FRACTION 

OF  THE 

FIRST 

YEAR. 

c 
c 
c 
c 
c 

IF  (X.LT.1.0)  THEN 
CAP(I)=FCAP(1) 
IF  (METHOD. EQ.1.0R.X.GE.TREPL)  THEN 

DTF(I)=FDT(5) 
ELSE 

IF  (X.GE.TREPL)  THEN 

DTF(I)=FDT(1) 
ELSE 

DTF(I)=FDT(1)/2 
END  IF 
END  IF 
ELSE 

R=M-IT 

IF  (R.EQ.O.O)  M=M-1 
CAP(I)=FCAP(M+1) 
DTF(I)=(FDT(M)+FDT(M+1 ) )/2 
END  IF 
END  IF 
20  CONTINUE 
RETURN 
END 
C 
C 
C 
C 
C  *»*«*»»»»»»»*»»»»» 

C  »SUBROUTINE  DTIME* 


SUBROUTINE  DTIME  (Y,DTY) 
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COMMON  DTF( 1000), CAP( 1000), ENA,ENB,OMA, OMB, REV, SRC, FOR, TREV,TLOOP, 
1T0TALA,T0TALB,SAV,BCR,T,FCR, LAMBDA, TREPL,QEN,QOM,QCOMP,QSSV,QTEMP, 
1QRM,QST0R,QFSV,QPCSV,QSPR,  DELTA,  SPRMW.MW.TI,  STEP,  FRAC,STOR,COMP,SS 
1V,TTEMP,RMA,RMB,FSV,PCSV,SPR,NI, NIG, IT, METHOD, MODEL, NP,N,ICOUNT,TE 
1MP 

REAL  LAMBDA 
TEX=(TREPL-TTEMP)/2 
C1=TTEMP/TREPL 
C2=1-C1 

C3=( 1/( 1+LAMBDA*TREPL) ) **(N-1 ) 
IF   (METH0D.EQ.1)    THEN 

DTY=C2+(1-C3)*C1 
ELSE 

X=Y+TEX 

CALL  GOT  (X.GTX) 
AA=1-(GTX«»(N-D) 
DTY=C2+AA*C1 
END  IF 
RETURN 
END 
C 
C 
C 
C 

Q  **•••* 

C  «DATA* 

SENTRY 

LAMBDA, REPAIR/REPLACE, TEMPORARY   REPLACEMENT, SPARE  MW  CAPACITY 

0.010     0.739726     0.616438       500.000 
//PLANTS, #INTERVALS(EXCEPT  GOT)  ,#INTERVALS(G0T)  ,METH0D,M0DEL 

13        10        10  1  4 

SPARE  COST, SPARE  SALVAGE  VALUE, SPARE  STORAGE 

2200.000      1000.000  20.000 

MAXIMUM  EVALUATION  TIME(MET, INTEGER) ,MET(REAL) , TIME  STEP, POOL  LIFETIME 

15  15.000  15.000  15.000 

ENERGY, O&M, NEW  COMPOSED  COMP  SAL, FAILED  COMP  SAL, TEMPORARY, RESERVE  MAR. 

74.898       200.000     2770.000       400.000  40.000  40.000        186.000 

PLANT  MW  CAPACITY, FCR, FOR, SRC, REV(YRS) 

800.000  0.174  0.220      1000.000  25.000 
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ABSTRACT 

The  purpose  of  this  study  was  to  refine  techniques  developed  in  a 
previous  study  on  the  economic  feasibility  of  sharing  a  spare  critical 
component  among  several  electric  generating  stations.   Total  costs 
attributable  to  a  particular  plant  component  were  estimated  under  two 
component-management  plans :   Plan  A  in  which  no  preparations  are  made 
prior  to  a  component  failure,  and  Plan  B  in  which  a  spare  component  is 
purchased  and  made  available  for  temporary  use  should  a  plant  component 
fail.   The  modeling  techniques  previously  developed  to  estimate  the 
difference  between  Plan  A  and  Plan  B  costs  over  the  lifetime  of  the 
spare-component  pool  were  studied  and  refined  in  two  major  areas.   The 
first  refinement  was  to  the  failure  model  that  describes  both  the 
probability  of  spare  availability  and  the  fraction  of  the  component 
repair/replacement  time  that  the  spare  is  not  used  as  a  temporary 
substitute.   A  Monte  Carlo  simulation  code  was  used  to  investigate  three 
different  stochastic  models.   A  fourth  failure  model,  based  on  a  renewal 
theory  description,  was  also  studied.  The  second  area  of  refinement  was 
to  the  economic  model  so  that  costs  could  be  estimated  on  a  continuous- 
time  basis  rather  than  at  year  ends  as  was  done  in  the  earlier  study. 

A  computer  code  was  written  that  incorporated  the  continuous-time 
economic  models  and  the  four  different  stochastic  failure  models.   This 
code  was  used  to  investigate  the  sensitivity  of  various  assumptions 
about  the  spare-component  pool  to  the  benefits  afforded  by  the  pool. 
Results  from  the  program  were  also  compared  to  results  obtained  in  the 
earlier  study.   From  this  investigation,  it  was  found  that  for  a  system 
of  five  or  more  plants,  consideration  of  the  possibility  of  the  spare 


failing  while  used  as  a  temporary  substitute  had  little  effect  on  the 
benefits  of  a  spare  component  pool.   Of  greater  importance  was  the 
choice  of  model  to  describe  the  spare  availability  and  the  duration  that 
the  spare  is  to  be  used  as  a  temporary  substitute.   For  long  spare  pool 
lifetimes,  constant  asymptotic  failure  models  gave  excellent  results; 
however,  for  short  pool  lifetimes  it  is  important  to  use  a  detailed 
description  of  spare  availability. 


